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FIGURE  la.  System  detection  probability  for  one  area  of  regard  (J=  1 )  of  rough  aluminum  as  a 
function  of  number  of  measurements  (AO,  prior  probability  P(H\),  and  the  pseudo  SNR  parameter 
y .  Contour  surfaces  show  where  Pdetect(system)  =  constant.  Three  contour  surfaces  are  shown  at 
Pdetect( system)  values  of  0.2  (blue),  0.5  (green/yellow),  and  0.8  (red) . 20 

FIGURE  lb.  Pseudo-STVT?  parameter  ;/  as  a  function  of  the  mean  fill  factor,  E(f),  and  the  mean 
contamination  thickness,  E(h),  where  the  surface  contamination  density  is  G(jug / cm2)  =  234  x  h 


where  h  is  in  microns) . 20 

FIGURE  2a.  An  expanded  view  of  Fig.  la  to  show  the  behavior  for  y< 4 . 21 


FIGURE  2b.  An  expanded  view  of  Fig.  lb  to  show  details  of  pseudo  SNR  parameter  y  as  a 
function  of  fill  factor  and  contamination  thickness.  Combinations  of  E (h)  and  E (/)  resulting  in 


y> 4  are  shown  as  dark  red  pixels,  y< 4  are  mapped  to  other  colors  as  shown . 21 

FIGURE  3.  Same  as  Fig.  la  but  J=3  areas  of  regard  are  used  for  detection  in  (19) . 22 
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FIGURE  4.  Detection  probability  of  local  detectors  Pdetect  .  i  n  (19)  for  simulation  1  (footprint 

contamination  Fig.  Al)  as  a  function  of  the  2nd  binary  threshold rj2] . .  j=\  is  AOR  for  rough 

aluminum,  and  j= 2  is  AOR  for  painted  aluminum  surface.  P(//\)=0.4.  The  system  detection 
probability,  under  the  condition  that  ^21  =  r/22 ,  is  also  shown  in  black.  Pdetecij=2  reaches  a 

maximum  of  0.926  at  r/22  =  2 .  Red  dot  shows  the  largest  r/2l  such  that  PMM  j=l  exceeds  a 
detection  probability  of  0.999  (occurs  at  rj2l  =3).  Black  dot  shows  the  location  where  the  system 
detection  probability  exceeds  0.999  (occurs  at  ij2X  =r/22  =5).  Combining  the  two  regions  allows 
larger  r/2J  thresholds  to  achieve  the  same  detection  probability  and  will  result  in  more  robust 
performance . 24 

FIGURE  5.  Overall  system  detection  probability  Pdetect(system)  in  (19)  for  case  1,  as  a  function  of 
all  possible  combinations  of  2nd  thresholds  rj2j Tn  (18).  P(H\)=0A.  rj2l  is  2nd  threshold  for  j=  1 
(rough  aluminum)  and  /;22  is  2nd  threshold  for  j= 2  (painted  aluminum).  Black  line  occurring  along 
the  diagonal  (where  =  >j22)  is  the  same  black  line  appearing  in  Fig.  6;  red  and  green  lines 
appearing  at  r/22  =  N  and  tj2l  =  N ,  respectively,  are  the  same  as  in  Fig.  4  (at  ij2j  =  N  the  detector 
will  never  alann  and  thus  information  from  the  jth  region  is  discarded).  Maximum  of  Pdetect(system) 
occurs  at  r)2X  =  1 ,  rj22  =  2  at  a  value  of  0.9999996  . 24 

FIGURE  6.  Same  as  Fig  4  but  for  contamination  footprint  (case  2)  of  Fig.  A2.  P(H\)=0.01 . 
Individual  detection  probability  PiaM  j  attain  a  maximum  of  0.934  and  0.919  for  j=l  and  j=2 

respectively  at  ri2j  =1 .  Maximum  detection  probability  for  the  system  is  0.995  at  ii2j  =1  . 25 

FIGURE  7.  Same  as  Fig.  5  but  for  contamination  footprint  (case  2)  of  Fig.  A2.  P(//i)=0.07  ....26 

FIGURE  Al.  Simulation  1  (case  1):  (left)  value  of/for  each  pixel  on  the  simulated  surface.  The 
sensor  FOV  (blue  box,  1/1 0th  the  linear  dimension  of  the  image)  moves  over  the  contaminated 
surface,  (right)  The  fill  factor  as  seen  by  the  FOV  for  the  contaminated  surface.  The  likelihood 
for  the  sensor  FOV  to  encounter  contamination  is  P(Hi)=0.4 . 30 

FIGURE  A2.  Simulation  2  (case  2):  same  as  Fig.  Al  except  that  individual  contamination 
events  are  larger  in  cross-sectional  area;  it  takes  fewer  events  for  the  same  total  area  of 
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Abstract 

We  developed  a  framework  to  evaluate  the  performance  of  a  detection  system  for  contaminated 
surfaces.  We  employ  the  radiative  transfer  model  for  contaminated  surfaces  (Ben-David  and 
Davidson,  ECBC-TR-1084,  2013)  and  transform  the  physical  model  into  a  stochastic  probability 
model  with  which  detection  probability  and  false  alarms  can  be  estimated  for  scenarios  of 
interest.  Our  algorithm  employs  a  data  fusion  approach  known  as  a  distributed  binary  integration 
system  (also  known  as  a  double-threshold  detector,  or  m-out-of-n  detector)  in  order  to  combine 
the  individual  detection  results  from  multiple  scans  over  several  potentially  contaminated  areas. 
With  our  probability  model  we  can  explore  the  parameter  space  (e.g.,  number  of  measurements, 
time  to  detect,  area  to  monitor,  sparsity  of  the  contamination,  field  of  view,  etc.)  and  study  the 
tradeoffs  between  parameters  that  affect  the  overall  system  detection  performance.  We  can  also 
use  the  stochastic  model  to  set  sensor  requirements  for  a  contamination  scenario.  We  presented 
plots  that  demonstrate  the  interaction  between  parameters  and  an  example  for  the  detection  of  a 
potassium  chlorate  contaminated  “car”  with  a  CCfi  tunable  laser  system. 

Subject  Terms  Radiative  transfer,  contaminated  surfaces,  detection,  rough  surface,  BRDF, 
reflectance,  fill  factor,  distributed  binary  integration,  CFAR,  data  fusion,  distributed  sensor 
system,  double-threshold  detector,  m-out-of-n  detector,  potassium  chlorate,  probability  theory, 
system  performance,  probability  of  detection  and  false  alarm. 

1.  Introduction 

In  a  recent  publication1  we  have  developed  a  radiative  transfer  model  with  which  we 
addressed  long-wave  infrared  (LWIR)  passive  and  active  spectral  reflectance  measurements  of 
potassium  chlorate  or  ammonium  nitrate  contaminated  rough  surfaces.  The  radiative  transfer 
model  was  developed  from  physics-based  principles  with  the  aid  of  empirical  approximation  of 
the  uncontaminated  rough  surface  reflectance.  In  this  report  our  objective  is  to  develop  a 
framework  for  the  radiative  transfer  model  where  the  parameters  are  allowed  to  be  stochastic 
(i.e.,  fluctuate)  and  thus  the  new  model  will  be  able  to  predict  the  performance  of  a  notional 
system  with  traditional  figures  of  merit:  probabilities  of  detection  and  false  alarm.  In  this  report 
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we  only  address  an  active  infrared  laser  system  (pulsed  or  continuous  wave)  that  probes  the 
contaminated  surface  and  measures  the  reflectance  at  the  specular  reflectance  angle  (e.g.,  a  bi¬ 
static  lidar  system)  or  at  the  backscatter  angle.  In  future  reports  we  will  enlarge  the  framework  to 
passive  remote  sensing  where  the  reflected  thermal  radiation  and  for  a  Raman  system  that  is  used 
to  interrogate  contaminated  rough  surfaces. 

The  main  steps  in  developing  a  framework  for  stochastic  radiative  transfer  model  are: 

I.  Identifying  the  relevant  physical  parameters  to  become  stochastic  random  variables  (RV). 

II.  Choosing  probability  models  for  the  selected  RVs. 

III.  Combining  the  RVs  with  the  physical  radiative  transfer  model  into  a  probability  model 
for  the  measured  (predicted)  signal. 

IV.  Selecting  a  signal  processing  method  with  which  spectral  signal  is  transformed  to  a 
detection  score. 

V.  Applying  data  acquisition  constraints  to  the  detection  scenario:  time  to  detect,  surface 
area  to  cover. 

VI.  Determine  system-level  performance:  combining  many  single  (local)  measurement 

probabilities  of  detection  and  false  alarm  into  a  scenario  (global)  probability  of  detection 
and  false  alann. 

2.  Overview 

We  start  by  giving  a  brief  overview  of  a  practical  scenario  for  which  we  envision  the 
implementation  of  our  stochastic  model.  A  suspected  car  at  a  check-point  is  scanned  with  a 
LWIR  laser  beam  and  reflected  signals  are  measured  throughout  key  areas  of  the  suspected  car. 
Within  a  given  time  and  after  many  measurements  are  acquired-where  each  measurement  is  a 
vector  with p  spectral  bands-a  detection  score  is  broadcasted  by  the  algorithm.  The  detection 
score  is  evaluated  with  respect  to  a  probability  of  false  alann  (the  car  is  falsely  declared  as 
“contaminated”)  and  the  probability  of  missed  detection  (equal  to  one  minus  the  probability  of 
detection)  where  the  car  is  contaminated  but  the  detection  algorithm  did  not  broadcast  the  proper 
“contamination  detected”  declaration.  The  parameters  “missed  detection”  and  “false  alann”  drive 
the  detection  algorithm  in  opposite  directions;  when  the  former  decreases  the  later  increases.  In 
the  scenario  of  a  suspected  car  we  are  willing  to  endure  higher  false  alarm  rates  in  order  to 
minimize  the  chance  of  a  missed  detection,  because  we  assume  false  alarms  will  only  necessitate 
a  secondary  (e.g.,  manual)  investigation  of  the  car  at  the  check-point  location. 

The  basic  radiative  transfer  model  [Eq.  10,  in  reference  1]  for  the  reflected  signal  from  a 
contaminated  surface,  illuminated  at  incidence  angle  9t  and  measured  at  reflectance  angle  6  is 
given  by 


2 


R{9,9i)  =  {\-f)xR0+fx 


-h\ 


Rt+(l-R,)R0xe 


1  1 
- + - 

cos(#,)  cos(0) 


(a+gAn) 


R0  =  p0  exp(a0T  2  +  b0A  1  +  c0) 
R,  =  Aexp(a,r2  +b,Xi  +c,) 


(1) 


where, /is  the  fill  factor  of  the  contaminated  area  (the  fraction  of  illuminated  spot  that  is 
contaminated),  h  is  the  contamination  thickness  (e.g.,  modeled  as  a  film),  Rl)(0,0i)  is  rough 

surface  reflectance  from  the  car,  Rt{ 0,0 /  is  the  reflectance  of  the  contamination  (i.e.,  target 
material)  modeled  as  a  rough  surface  due  to  morphology  of  the  contamination,  /  and  9  are  the 
incident  and  reflected  angles  measured  counter-clockwise  from  the  surface’s  normal  (i.e.,  for 
specular  reflectance  geometry  9  =  -9i ,  and  for  lidar  backscatter  geometry  0  =  9 \ ),  a  is  the 
spectral  absorption  coefficient,  and  gAn  is  a  dispersion  correction  to  a  .  In  (1),  top  line,  R0 ,  Rt, 
and  R  are  wavelength  dependent  (wavelength  dependence,!  ,  is  omitted  for  clarity);/  h,  a,  b,  c, 
are  all  scalars.  Both  Ro  and  Rt  are  given  with  empirical  coefficients  {a,  b,  c )  and  the  spectral 
Fresnel  reflectance  coefficients  p0  and  pt  are  computed  at  incidence  angle  0; .  In  [  1  ]  we  showed 

very  good  results  of  applying  the  model  (2)  for  specular  geometry.  In  this  report  we  speculate 
that  also  for  lidar  backscatter  geometry  (that  can  be  implemented  with  a  hand-held  sensor)  the 
model  given  by  (2)  provides  good  results. 

Thus,  given  atmosphere  with  transmission  talm  {range)  where  range  is  the  folded  path 
from  the  source  to  the  surface  and  back  to  the  detector  an  illuminating  source  with  constant 
strength  Lsource{0i )  at  a  distance  from  the  surface  (distance  rsource )  and  a  detector  at  a  distance 

Selector  anc'  angle  0  from  the  surface,  the  measured  reflected  signal  is  given  by 

M, „((?)  =  !„„,(«,)- - (= - vR(ff„g)  =  K„xR(g,ffl)  (2) 

where  Ksystem  =  Lsource{9t)talm{rsource  +  rdetector)“2  is  a  system  constant  (for  a  given  geometry)  and 
R(0,0t)  is  given  in  Eq.  (1). 

To  develop  a  probability  model  for  Mdetector  we  need  to  develop  (or  assume)  a  statistical 
model  for  all  random  variables  (RV)  in  the  model  equation  and  then  modify  the  statistics  of 
R{0,0t)  with  the  scaling  factor  K  tem  .  If  needed,  one  can  regard  Ksystem  as  a  RV  (e.g.,  due  to 

fluctuations  in  atmospheric  transmission  or  laser  power)  -  in  this  work  (for  simplicity)  K  tem  is 
not  a  RV.  The  objective  is  to  develop  a  probability  density  function  {pdf)  for  the  at-aperture 
signal  Mdetector,  and  use  it  in  a  detection  algorithm  (e.g.,  a  matched  filter  detector)  to  produce 
detection  scores  that  can  be  evaluated  for  probability  of  detection,  probability  of  false  alarm,  and 
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probability  of  missed  detection.  In  finding  probability  model  for  the  RVs  in  (2)  we  strive  for  a 
balance  between  complexity  of  the  statistical  description  and  ease  of  use  at  the  expense  of 
accuracy.  It  is  counterproductive  to  describe  RVs  with  complex  mathematical  functions  that  we 
cannot  use  when  we  need  to  manipulate  RVs  in  the  form  of  multiplications  and  additions  of  RVs 
that  lead  to  the  pdf  of  (2). 

The  RVs  that  we  address  in  this  study  are  the  fill  factor/  contamination  thickness  h, 
rough  surface  reflectance  Rq,  and  the  reflectance  of  the  contaminated  target,  Rt.  The  statistics  of 
the  fill  factor/is  affected  by  the  field  of  view  (FOV),  the  cross-sectional  area  distribution  (pa{a) 
of  the  contamination,  and  the  total  mass  of  the  contamination  per  unit  area.  The  spill  area- 
distribution  <pa(a)  is  the  footprint  distribution  of  the  contamination  on  the  surface.  While  (pa  ( a ) 
is  related  to  the  intrinsic  diameter  distribution  function  <f>  diameter  ( diameter)  of  the  contamination 
particles,  the  two  are  not  the  same.  The  spill-distribution  (pa  (a)  describes  the  interaction  between 
single  particles  due  to  conglomeration,  clumping,  how  particles  fall  out  (spill)  of  a  sack  that 
moved  by  a  car,  etc.,  whereas  ®  diameler(diameter)  is  the  size  distribution  of  the  individual 
particles  before  they  were  spread  (spilled)  on  the  surface  (e.g.,  a  size  distribution  of  a  powder 
before  it  was  put  in  a  sack).  The  contamination  footprint  cross-sectional  area  distribution  (pa  ( a ) 
is  used  to  compute  the  pdf  for  f  Pf  (/) .  The  particle  size  distribution  function 
O  diameter  (diameter)  may  be  used  to  compute  the  pdf  of  the  target  thickness  h,  Ph(h) ,  when  one 

assumes  a  single  layer  of  particles  on  the  surface  for  which  h  =  diameter .  In  this  work  we  regard 
the  pdfs  for / and  h  to  be  independent  of  each  other  (for  convenience  and  because  we  do  not 
know  the  mechanism  of  spillage  that  create/ and  h),  however,  if  one  choose  to  impose 
dependency  between /  and  h,  then  a  bivariate  Ph  f  (/?,  /)  can  be  used  (a  bivariate  pdf  for  two 

marginal  nonnals  is  easily  obtained,  and  for  non-nonnal  the  concept  of  copula  can  be  employed). 
In  the  notation  for  a  general  pdf  e.g.,  Pf  x  \  z ) ,  the  subscript  denotes  the  parameter  (y)  for  the 

pdf  and  in  parentheses  it  is  a  dummy  argument  (x).  The  notation  x  |  z  reads  as  x  conditional  to  a 
given  z  value;  Py(x)  denotes  an  unconditional  pdf. 

Let’s  assume  that  we  are  able  to  derive  the  pdf  Pr(R )  for  the  measured  reflectance  R  in 
(1);  then  with  a  simple  transfonnation  the  pdf  of  (2)  is  PM  (M)  =  K~ltemPR  (M  / Ksystem ) .  With 

PM  ( M )  we  can  formulate  our  binary  Neyman-Pearson  hypothesis  testing2  regarding  the  data  M 
under  Ho  and  ll\  scenarios,  where  for  the  null  hypothesis  Ho  the  interpretation  is  that  the  data 
does  not  contain  the  signal  of  interest  (e.g.,  a  target),  and  the  alternative  hypothesis  Hx  is  that  the 
data  does  contain  the  signal  of  interest  in  addition  to  noise  and  interferences.  Under  the  two 
hypotheses 

\M  I  H0  =M(fh  =  0)] 

\m  |  Hl  -  M(fh  ^  0)  J 
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The  fluctuations  in  M  are  due  to  both  additive  noise  (whose  source  can  be — but  not  limited  to — 
detector  noise),  and  to  multiplicative  noise  (e.g.,  due  to  clutter).  In  most  sensors  system  noise  is 
small,  and  most  of  the  noise  is  due  to  clutter  which  in  our  model  is  represented  by  the 
uncontaminated  rough  surface  pdf  PK  ( R0 ) . 

Analyzing  detection  performance  requires  selecting  a  signal  processing  method  with 
which  to  determine  the  presence  of  the  target  signals  in  the  data.  We  have  chosen  the  standard 
matched  filter  algorithm,  which  is  optimal  for  additive  Gaussian  noise.  Although  our  problem  is 
not  described  by  a  simple  additive  model  (and  thus  the  matched  filter  is  not  optimal),  the 
matched  filter  is  a  standard  technique  used  as  a  benchmark  for  comparing  algorithms,  leads  to 
convenient  performance  analysis,  and  is  therefore  a  suitable  choice.  The  matched  filter  detection 
score  is  given  by 

,core  [M  -  E(M  |  Hn)]T  cov(M  |  H0)~'  8 
Stcov(M\H0Y'S  y 

where  E(-)  is  a  mean  (expectation),  cov(-)  is  a  covariance  operation,  superscript  “7”  is  a  matrix 
transpose  operation  (when  M  is  a  column  spectral  vector  of  measurements  with  p  spectral  bands 
and  Mt  is  a  row  vector),  8  is  the  target  spectral  vector  that  represent  the  spectral  direction  of 
the  intrinsic  contamination  (to  be  discussed  in  section  5).  The  simplest,  but  not  necessarily  best, 
choice  for  8  is  the  spectral  absorption  coefficient  er  in  (1).  For  the  detection  scores  we  compute  a 
probability  of  detection  Pdetect  and  probability  of  false  alarm  Pfals^alarm .  The  scores  in  (4)  are 

function  of  the  /!<>  condition,  namely,  a  function  of  the  local  rough  surface  properties.  Thus,  it  is 
prudent  to  assume  that  different  parts  of  the  suspected  car  will  have  different  H0  signals  (e.g., 
due  to  non-uniformity  of  materials,  and  geometry  of  orientation  with  respect  to  incidence  angle). 
In  the  general  detection  case  we  have  JAOR  areas  of  regard  (AOR)  in  the  suspected  car,  each  area 

of  regard  produces  Nscores.  The  scores  for  each  J AOR  are  characterized  by  Pdetect  and  Pfalse_alarm . 

We  combine  the  JAOR  x  Nscores  local  scores  and  probabilities  into  one  (global)  set  of  probabilities, 
^detect (car)  and  PUse_ahvm(car)  with  the  “binary  integration”  method3'4,  a  method  that  is  used  in 
radar  detection  algorithms.  The  number  of  AORs  can  be  set  according  to  distinct  areas  that  affect 
the  substrate  reflectivity  R  |  H0 :  the  number  of  different  materials  of  the  uncontaminated  car’s 

surfaces  (e.g.,  rubber,  aluminum,  painted  aluminum,  etc.),  and/or  the  orientation  of  the  surfaces 
(e.g.,  horizontal  hood  and  trunk,  vertical  doors),  and  could  also  depend  on  a  priori  information 
regarding  the  likelihood  of  different  areas  of  the  car  to  be  contaminated. 

The  total  number  of  detection  scores  J AOR  x  Nscores  is  a  function  of  external  constraints 
such  as  total  time  ( T)  available  for  scanning  a  suspected  car,  sensor  acquisition  rate  (laser  pulse 
repetition  frequency,  PRF,  for  an  active  sensor),  and  the  sensor  field  of  view  (solid  angle,  FOV 
in  steradians)  whose  viewing  area  is  r^ourceFOV  .  For  convenience,  we  assume  that  the  laser  beam 

divergence  matches  the  field  of  view,  hence,  rjourceFOV  is  the  illuminated  spot-area  on  the 
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contaminated  surface.  The  internal  relationship  between  the  external  constraints  and  a  simple 
flowchart  is  given  by 

external  ]  ,  .  ,  [detection 

constraints  II.  .1  iscores 


J  AOR  X  N scores  =  T  X  PRF 


N.. 


A/J 
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r  FOV 

source 


detector 
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>  < 
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>  < 
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detection  O AOR  ) 
false-alarm  O  AOR  ) 


overall 

probability 
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false-alarm 


(5) 


J 

Note  that  in  (5)  A  is  the  sum  of  the  viewed  areas  for  each  JAOR  areas  of  regard,  A  =  I A  • 

3= 1 

Nscores  is  the  number  of  measurements  for  each  area  of  regard  that  is  collected  with  a  sensor  with 


a  given  FOV (sr),  and  that p  spectral  bands  are  collected  every  PRF  1  seconds.  The  illuminated 
area  (by  the  active  source)  for  a  single  laser  shot  on  the  contaminated  surface  is  rs2ourceFOV  .  Solid 
angle  FOV  (sr)  is  the  square  of  the  linear  FOV  in  radians.  In  (5)  each  of  the  AORs 
contains  Nscores,  however,  one  can  allow  a  different  number  of  measurements  for  each  j,  Nj ,  and 

j 

that  JAOR  x  Nscores  =  ^  J ;  /V j  .  The  sensor  can  be  a  single  pixel  detector  or  an  array  for  which  the 

j=i 

acquisition  of  Nscores  may  take  less  time.  For  example,  for  N=\  00,  J=  6,  and  7=3 Os,  the 


PRF  =  J  xN /T  =  20 Hz  ,  which  implies  that  a  vector  with p  spectral  bands  must  be  acquired 
every  5 ms  -  a  difficult  task  not  to  be  underestimated.  With  the  snapshot-advantage  architecture5 
(see  Fig.  1  in  reference  5)  this  task  can  be  achieved;  one  can  have  a  detector  with  x-by-v-by-k 
data  cube  where  x  =  y  =  V77 ,  and  k  =p,  that  is  acquired  simultaneously  for  each  of  the  J  areas  of 
regard.  In  order  for  an  active  sensor  to  take  advantage  of  the  snapshot-advantage  architecture 
(i.e.,  an  imaging  receiver)  one  would  need  multiple  sources  or  else  ensure  that  the  illuminated 
spot  size  is  large  enough  to  illuminate  each  pixel  in  the  image  simultaneously.  In  this  report  we 
assume  that  N  and  J  are  the  external  parameters  to  the  global  probability  model  and  we  gloss 
over  the  acquisition  strategy  that  was  employed  to  obtain  them. 


3.  Probability  model  for  stochastic  variables 

We  describe  the  options  for  the  principal  RVs  and  our  choice  for  a  statistical  model.  When 
selecting  a  statistical  model  for  a  RV  we  consider  the  range  of  the  RV:  e.g.,  the  support  range  for 
contamination  thickness,  h,  is  zero  to  infinity  or  zero  to  hm ax;  the  support  range  for  a  fill  factor,/ 
is  zero  to  one.  There  is  more  than  one  pdf  that  can  describe  the  range  of  values  of  a  given  RV. 

The  most  likely  choice  for  a  pdf  in  the  absence  of  any  prior  knowledge,  is  a  pdf  whose  entropy  is 
a  maximum.  For  example,  given  a  RV  x  that  is  known  to  exist  only  within  a  limited  range 
(support)  Xj  <  x  <  x2  and  is  zero  outside  the  interval  [Xj  ,x2  ]  the  most  likely  pdf  is  a  uniform  pdf 
i.e.,  with  our  maximum  ignorance  we  must  accept  that  all  values  x,  <  x  <  x2  are  equally  likely.  If 
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we  have  prior  knowledge  about  x  that  the  mean  E(x)>  0 ,  and  the  support  0  <  x  <  oo  ,  then  an 
exponential  pdf  is  the  most  likely  pdf  for  x;  If  we  know  the  mean  E(x)  and  variance  V  (x)  of  x, 
then  the  nonnal  distribution  for  which  -  oo  <  x  <  oo  is  the  most  likely  pdf 

In  this  study  we  assume  measurements  are  collected  at  p  wavelengths,  hence,  Ro,  Rt,  and 
Mare  multivariate  RVs  (i.e.,/?-by-l  vectors),  whereas / and  h  are  univariate  RVs  (i.e.,  1-by-l 
scalars). 

3.1  Surface  roughness  Ro  and  Rt 

The  RV  for  surface  roughness  R0  =  pQ  exp(a0/T2  +  b0A~'  +  c(j )  =  p0D  where 
D  =  exp(a0/T2  +  b0A  1  +  c(} )  is  the  decay  function  (see  equation  3  in  reference  1)  contains  three 
RVs,  -oo  <  (a0,b0,c0)  <  oo  and  is  bounded  in  the  support  range  0  <  R0  <  1 .  If  we  knew  the 
means  (pa0 ,  pb0 ,  pc0 )  and  the  variances  (<r20 ,  af ,  <r20 )  of  (a0 ,  b0 ,  c0 ) ;  we  can  assume  a  nonnal 
distribution,  N(p,a2 ) ,  as  the  most  likely  pdf  for  (a0,b0,c0),  and  R0  is  given  by  lognonnal  pdf 
LN(p,a 2)  as 

Rq  ~  LN(\n(p0 )  +  pD ,  2/j) 

PD  —  RaO^-  +  RbO^  RcO 

<  2,  =  af(r2xr2)T  +  cr^XA-y  +  a20llT 
D~LN{pDXD ) 

where  1  is  a  p- by-1  vector  with  all  ones.  Note,  that  in  (6)  we  did  not  enforce  the  upper  bound 
constraint  R0  <  1  (the  lower  bound  is  assured  to  be  greater  than  zero  when  using  lognormal  pdf) 

Nevertheless,  since  we  know  from  lab  experiments  that  D ~  LN(p/),YJf))  is  spectrally  smooth 
function  due  to  the  spectral  dependence  imposed  by  the  2nd  order  polynomial  a„A  2  +  h,f  1  +  c0 , 
we  can  choose  scalar  values  for  ( a0,b0,c0 )  to  obtain  a  reasonable  decay  function!),  and 
furthermore  we  can  ensure  that  our  choice  for  ( pD,HD )  will  produce  a  physical  rough  surface 
reflectance  R0<  1 .  The  same  can  be  done  to  produce  a  pdf  for  the  contamination  (viewed  as  a 
rough  surface)  R,  using  (a, ,bt,ct ) , 

Rt  ~  LN (ln(p, )  +  pDt ,  X  !)t ) 

Mot  =  RaX~  +  RbX  +  Met 

<  =  °IM-2XI-2)T  +  ai() +  <rl  11T 

D,  ~  LN(pDt,ZDt) 

In  future  work  we  may  generalize  (6,  7)  and  treat  the  Fresnel  reflectance  p0  and  pt  as  RVs, 
which  would  model  a  change  in  composition  of  the  car’s  surface  and  the  target. 
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3.2  Contamination  thickness  h 

In  our  physical  model  the  thickness  h  of  the  contamination  is  an  average  value  across  the 
contaminated  areas  within  the  field  of  view.  Due  to  the  non-uniformity  of  contamination  we 
expect  variations  in  h,  hence  fluctuations  in  the  transmission  t  x  exp  {-ah)  through  the 
contamination.  Our  first  option  is  a  uniform  pdf  for  h.  We  can  assume  that  in  practice  h  is 
bounded,  hmin  <h<  /?max ,  where  himx  is  chosen  such  that  the  transmission  t  x  exp  (-ah)  is 
practically  zero  (i.e.,  within  sensor  noise).  A  uniform  pdf  for  0<  h<  hmiiX  has  maximum  entropy 

and  thus  is  advantageous  (maximum  entropy  implies  high  likelihood  of  occurance).  However, 
the  pdf  of  transmission  (easily  obtained  with  transformation  of  variables  for  the  pdf  of  h)  will 
have  higher  occurrence  for  low  transmission  than  high  transmission,  a  scenario  that  we  do  not 
think  as  being  reasonable.  On  the  contrary,  we  think  that  most  transmission  events  across  the 
contaminated  surface  will  be  of  high  value  (i.e.,  associated  with  small  h).  Our  second  choice  is 
an  exponential  pdf  for  h.  We  know  that  h>  0  and  thus  (for  maximum  entropy)  we  are  tempted  to 
choose  an  exponential  pdf  for  0  <  h  <  oo  ,  hence  larger  probability  for  small  contamination  (thin 
layer  h)  than  for  large  amount  of  contamination.  The  problem  with  exponential  pdf  for  h  is  its 
high  non-zero  density  for  h= 0  and  thus  the  exponential  choice  imposes  restriction  on  the 
interpretation  of  the  H\  scenario  (contamination  is  present),  i.e.,  Ph(h=0\Hi)  =  0  is  required  by 
definition  of  H\.  Our  third  option  is  a  simple  linear  pdf  pdf  (h)  x  const  /  h  for  A  between 


hmin  <h<  /zmax  such  that  pdf[hmin )  >  pdffimax )  (an  additional  infonnation  to  set  the  slope  of  the 

pdf  is  required)  and  hence  the  pdf  of  the  transmission  will  have  higher  occurrence  for  small  h 
than  for  large  h  (the  behavior  we  desire).  Other  options  for  pdfs  are  the  use  of  complex  (but 
flexible)  functions  such  as  beta  and  gamma  pdfs  that  can  be  constrained  (with  additional 
information)  to  give  the  expected  behavior  for  the  pdfs  of  h  and  the  transmission.  If  we  have 
information  about  the  size  distribution  ® diameter  {diameter)  of  the  particles,  we  can  use  it  to  form 

a  pdf  for  h.  For  example,  we  can  assume  that  the  contaminated  surface  is  composed  of  a  single 
layer  of  particles,  and  thus,  the  thickness  is  proportional  to  a  diameter.  Thus,  the  pdf  of  h  will 
have  a  similar  shape  as  the  particle  size  distribution.  The  experimental  data6  suggests  that  a 
monolayer  of  contamination  is  reasonable.  Our,  choice  for  pdf  for  h  is  to  use  the  size  distribution 
diameter  (diameter)  or  if  is  not  available  to  use  the  simple  linear  pdf  (li)  x  const  /  h  ,  for  h 


between  hmm  <h<  hnrdX ,  that  is,  Pfh) 
density  at  /zmax. 


^max  ~  tl 
(^max  _  ^min  ) 


which  has  a  mode  at  hm ;n  and  zero 


3.3  Contamination  fill  factor  f 

A  pdf  for  /  was  not  given  in  the  MIT-LL  report6;  only  the  footprint  cross-sectional  area  of  the 
spill  distribution  cpfa )  was  given.  Thus,  we  had  to  compute  the  pdf  for  the  fill  factor.  We 
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compute  Pf(j)  with  two  methods:  via  simulations  and  with  a  theory  (see  Appendix  A  for 

details  of  the  simulation  method).  A  brief  description  of  the  two  methods  is  as  follows. 

(1)  a  Pj  (f)  is  computed  with  a  simulation  program  where  we  construct  a  random  surface 

populated  with  random  cross-sectional  footprints  that  are  sampled  from  (pa  (a)  and  are 

unifonnly  distributed  on  a  surface.  The  footprints  are  not  allowed  to  intersect  (i.e.,  no 
overlap  between  footprints).  With  this  complex  simulation  we  compute  a  fill  factor  for 
a  given  field  of  view  by  moving  the  FOV  over  the  simulated  random  surface  (a  digital 
image)  and  computing  the  area  of  all  the  particles  in  the  simulations.  Due  to  pixilation 
of  the  objects  in  the  digital  image  (the  size  of  FOV  and  the  spill  footprint)  the 
computed/)  (/)  may  contain  sharp  spikes;  therefore,  we  typically  smooth  the  pdf. 

(2)  We  compute  a  Pf(f)  with  a  theoretical  method  where  the  location  of  the  contamination 

sites  are  located  at  random  locations7  based  on  random  theory  of  distribution  of  random 
centers  in  a  hyper  volume  (for  our  application  we  reduced  the  volume  to  a  surface).  The 
spill  footprint  sizes  are  sampled  from  (pa  (a)  .  The  probability  of  intersection  between  a 

given  footprint  (modeled  as  a  circle)  and  the  edge  of  the  FOV  (modeled  as  a  rectangle) 
is  computed.  The  theoretical  method  is  based  on  detailed  geometrical  calculations  for 
intersection  of  shapes  ( FOV  and  footprint);  it  utilizes  approximations  that  are  only  valid 
under  certain  conditions. 

At  present,  we  choose  to  use  the  simulation  program  to  produce  fill  factor  probability.  In  future 
work  we  intend  to  improve  the  accuracy  of  the  theoretically  derived  pdf. 

The  prior  probability  P{HX)  for  II\  event  (i.e.,  the  likelihood  that  the  contamination  will 
be  present  within  an  FOV)  can  be  easily  computed  in  the  simulation  by  the  ratio  of  number  of 
FOVs  that  are  not  empty  (i.e.,  FOV s  that  contain  any  target  material)  to  the  number  of  FOVs  that 
are  contained  in  the  AOR.  For  example,  if  the  FOV  is  10  mr,  AOR  is  1  m  ,  and  the  distance  to 
the  AOR  is  10m,  then  the  number  of  independent  FOVs  is  1  m 2  /[1 0 mr  x  10] 2  =  100  .  If  in  the 
simulation  only,  say,  5  FOVs  contain  the  target,  then  P{HX)  =  5  / 100 .  The  prior  probability 
P(H\)  is  dependent  on  the  surface  density  of  the  contamination,  the  cross-sectional  area 
distribution  <pa(a) ,  and  the  size  of  the  FOV. 

4.  Statistics  of  the  signal  M 

It  is  extremely  difficult  to  derive  the  pdf  of  Mdetector  in  (2).  Our  approach  is  to  compute  the 
moments  of  Mdetector  and  then  to  fit  an  analytical  pdf  to  the  computed  moments.  Our  goal  is  to 
obtain  a  pdf  for  the  detection  scores  (4)  from  which  we  compute  a  probability  of  detection  and 
false  alarm.  Given  the  physical  model  (2)  Mdetector  =  Ksystem  x  R(0)  we  compute  its  nth  moment 

E{Mn)  =  KnsystemE[Rn{6)\  (8) 
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By  carefully  expanding  R"  (6 )  of  (2)  we  obtain  a  long  series  of  terms  that  is  an  «th  order 
polynomial  with  the  RVs  (f,  exp(-h),  R, ,  and  Rq)  where  each  term  of  the  series  is  of  the  form 
f"f  exp(-/j)"'‘  R"R'R"R0  where  the  sum  of  the  exponents  nf  +  nh  +  nRt  +  nR0  =  n  .  If  we  assume  that 

the  RVs  are  statistically  independent,  the  expectation  of  fn/  exp(-/?)"''  R"":  R'f"  is  simply  given 

by  exp(-/?)"'  R'f" )  =  E(f"f  )E(exp(-h)"h  )E(  Rf”  )E(  R'f" )  where  each  of  the 

expectations  can  be  computed  with  the  proper  statistical  model  chosen  for  that  specific  RV. 

With  first  and  second  moments  (n=  1  and  n= 2)  we  can  approximate  the  pdf  of  M  as 
multivariate  normal  (for  p  spectral  bands)  M  ~  N(jUm ,  ZM ) ,  or  a  lognormal  M  ~  LN(jUm  ,  ZM ) , 
where  the  population  parameters  are 

Mm=E(M ) 

ZM  =co v(M)  =  E(M2)-E(M)2  >  (9) 

=  E(MMt)  -  E(M)E(M)t 

The  advantage  of  choosing  a  lognormal  pdf  for  M  is  that  it  enforces  the  physical  constraint  that 
the  radiance  M  must  be  positive.  However,  normal  pdfs  are  easier  to  manipulate  mathematically 
(e.g.,  in  subtraction  operations  of  RVs)  than  lognormal  pdfs,  and  with  proper  choice  of  the 
population  parameters  for  normal  pdf  we  can  enforce  M  >  0  in  the  sense  that 

probability  pdf(M<0)  is  practically  zero  (due  to  the  fact  that  the  value  of  M  is  “far”  from  zero).  A 
more  complicated  pdf  (but  also  more  accurate)  for  M  can  be  given  with  the  Johnson  distribution 

o 

family  ,  by  a  fitting  procedure  that  uses  the  first  four  moments  of  M.  Our  choice  is  to  use  the 
simple  normal  pdf  for  M  for  which  the  statistics  of  the  detection  scores  (4)  are  more  easily 
obtained. 

5.  Statistics  of  the  detector  scores 

The  scores,  (4),  are  scalar  quantities  that  embody  a  spectrally- weighted  measure  of  how  similar 
the  measured  signal  is  to  the  expected  target  spectrum.  The  larger  the  score,  the  more  similar  the 
data  resembles  the  target  spectrum.  Given  the  detection  scores  it  is  of  interest  to  estimate  the 
probability  that  a  particular  score  value  can  be  attributed  to  an  Ho  scenario  or  to  an  Hi  scenario. 
Using  the  normal  approximation,  the  signal  Mis  distributed  by 

M0=M|i/0~V(//0,Z0) 

'  <u0  =  E(M  |  fh  =  0);  Z0  =  c ov(M  \fh=0) 

//,  =  E(M  |  fh  *  0);  Ej  =  co v(M  |  fh  *  0) 

where  (//„,//, ,  Z,P  Z, )  are  computed  from  the  first  and  second  moments  of  the  pdfs  of  the 
variables  (f,  h,  R0,  Rt).  The  target  direction  8  for  the  matched  filter  (see  (4))  is  chosen  to  be 
8  =  px-pa  where  we  note  that  //,  is  a  function  of  the  bare  surface  reflectivity  R0,  and  therefore 
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8  is  also  a  function  of  Rq.  In  practice,  /./,  may  not  be  known  a  priori,  hence,  one  can  use  a 


instead  of  8  -  a  choice  that  will  cause  a  drop  of  detection  perfonnance  that  can  be  estimated  with 
the  magnitude  of  the  cosine  angle  between  the  vector  X0°'5  a  and  the  vector  X0°'5(/h  -  Mo )  •  We 


note  that  the  drop  of  performance  affects  the  detection  probability  in  (12)  below,  which  is  the  1st 
stage  in  the  binary  integration  (to  be  discussed  in  section  6)  but,  can  be  somewhat  mitigated  at 
the  2nd  stage  of  the  binary  integration  where  the  sum  of  scores  above  a  given  threshold  is 
computed.  Equation  (4)  can  be  written  with  (10)  as 

[M-pJXf'8  1 


score  =  - 


sTxy  s 


S  =  jut -JU0 


•(H) 


If  Z0  is  known,  and  thus  is  not  treated  as  stochastic,  the  scores  are  distributed  with  nonnal 
statistics  as  two  normal  distributions  with  unequal  variances  located  at  zero  (for  H0 )  and  at  one 
(for  Hi)  given  by  ( 12) 


score\H0  ~  N(0,[ST  X0l  8]  l)  =  p0(x) 

score  |  Hx  ~  I,1  Xq1  S][St  X0'  Sf2)  =  P]  (x) 


r  = 


E {score  \Hl)~  E(score  \  H 0_)  =  ^0.5 


iJVar (score  \  H0) 

00 

PfioJM)  =  f  Po(x)dx  =  2-‘[l  -  erf(  J1  )] 

1  j2(8rZ018yI 


il -l 


^detection  (7)  =  f  Pi  (x)dx  =  2  ‘  [1  -  df  (  - - - - - - - 

i  y/2(8Tx0lx1x01s)(8TXo1dy 


yi 


(12) 


where  p0  and  p{  (top  two  lines  of  12)  are  the  pdfs  for  the  H0  scores  and  Hi  scores;  y  is  the 
distance  between  the  means  of  the  scores  under  Hi  and  Ho  in  units  of  standard  deviation  of  the 
Ho  scores  (note  that  for  X,  =  Z0 ,  /  =  SNR  of  the  H\  scores  and  is  a  sufficient  statistic  to 
completely  describe  detection  performance),  we  refer  to  y  as  “ pseudo  SNR”\  the  probabilities  of 
detection  PdctccUoli(p)  and  false  alann  /J;.ilse  ( q )  (bottom  two  lines)  for  a  given  threshold  //  are  the 

cumulative  density  functions  ( cdfs )  of  the  scores  and  are  given  with  the  error  function.  The 
parameter  y  is  related  to  the  separation  between  Ho  and  Hi  score  distributions,  but  if  X,  ^  X0  y 
does  not  completely  describe  the  separation  and  it  is  possible  for  two  scenarios  with  the  same/ 
exhibit  different  performance.  If  X0  is  not  known  a  priori  (e.g.,  we  are  asked  to  predict  the 
scores  distributions  at  a  later  day  for  which  we  do  not  have  Z„)  and  thus  is  a  random  variable, 
the  statistics  of  the  predicted  scores  are  given  by  confluent  hypergeometric  statistics9  when 
Xj  =x0;  when  X!  *  X0  the  distributions  in  reference  [9]  do  not  apply  and  must  be  modified. 
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6.  Combining  detector  scores  and  probabilities  from  multiple  scores  and  areas  of  regard 

In  the  general  scenario  we  interrogate  the  event  (a  suspicious  car)  in  J  areas  of  regard  (see 
external  constraints  in  (5))  where  for  the /h  area  we  have  A  detections  scores  with  probabilities 
of  false  alarm  and  detection  under  H0  and  Hi  (given  by  (1 1,  12)).  Each  of  the  J  AOR’s  can  be 
viewed  as  a  “local  detector”  that  collects  N  measurements  and  produces  N  scores.  The  challenge 
is  to  combine  all  the  N  scores  for  the  J  regions  ( J  x  N  results)  into  one  probability  of  detection 
(and  false  alarm)  for  the  event.  For  simplicity  and  ease  of  notation,  we  assume  an  equal  number 
of  measurements  AE  =  A  for  each  of  the  J  areas.  We  propose  to  use  the  radar  methodology3  of 

“binary  integration”,  also  known  as  “double-threshold  detector”,  and  “m-out-of-n  detector”. 
While  detectors  that  employ  binary  integration  techniques  (commonly  used  in  radar)  are  not 
optimal  in  perfonnance,  they  are  easy  to  implement  and  are  relatively  insensitive  to  a  single  large 
interference/clutter  that  might  exist  in  the  scene.  No  matter  what  is  the  energy  in  the  measured 
signal,  the  output  from  the  1st  threshold  is  a  binary  “1”  or  “0”,  and  thus  binary  integration 
detectors  are  robust  in  perfonnance  (less  susceptible  to  false  alarms)  when  the  background 
noise/clutter  is  a  non-Gaussian  statistics  with  high  tails  (reference  3,  p293-294)  that  can  produce 
large  scores  where  one  large  (false)  score  may  dominate  all  other  scores.  With  a  binary 
integration  method  all  scores  are  transformed  to  “1”  or  “0”  regardless  of  their  numerical  value.  In 
our  context,  the  binary  integration  addresses  the  problem  of  deciding  if  the  target  is  present  in  the 
/h  AOR  given  N.  Then,  we  address  the  following  question:  given  the  /  results  what  is  the 
probability  of  detection  for  the  global  event  (contaminated  car). 

From  a  system  perspective  (the  system  is  a  set  of  J  local  detectors,  each  of  which 
produces  N  measurements)  the  mathematical  statement  of  the  problem  is  as  follows.  For  the 
system  as  a  whole  the  objective  is  to  minimize  a  given  cost  function  for  the  global  event  and  to 
find  the  parameters  (the  two  thresholds)  for  each  of  the  /  detectors.  In  order  to  achieve  this 
objective  with  the  “binary  integration”  (double-threshold  detector)  method  we  need  to  set  two 
thresholds  for  each  of  the  /  local  detectors  (we  need  to  set  2  x  /  thresholds  for  the  system)  that 
together  will  satisfy  the  system  objective.  Fet  the  1st  threshold  for  the /h  local  detector  be  //, .  and 

the  2nd  threshold  be//2/ .  We  compare  the  zth  measurement  (out  of  N  measurements)  for  a  local 
detector  j  to  the  threshold  ^1;  and  produce  a  local  binary  score  ,s'/Y  =  1  when  score  M  >  ?jtj ,  or 
s ..  =  0  otherwise.  We  obtain  A  binary  scores  for  the /h  detector.  Then,  we  sum  the  A  binary 

N  N 

scores  and  compare  the  sum  ^  sp  to  the  2nd  threshold //2/ ,  Sf  =  YY.  >  //2/ ,  to  produce  an 

(=1  i=l 

integrated  detection  score  S  ■  that  is  one  when  the  sum  exceeds  the  threshold  value  (or  zero 

otherwise)  indicating  that  a  target  was  detected  by  jth  local  detector.  The  value  of  the  2nd 
threshold  72;  is  the  “/«”  in  the  name  “m-out-of-n  detector”  that  is  used  synonymously  with  the 

name  “binary  integration”.  The  sequence  of  A  binary  scores  (for  a  given  detector  j),  sJt ,  is  a 
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N 

Bernoulli  sequence,  and  the  sum  of  the  scores,  SJ  -  y  .y  n  is  a  binomial  random  variable.  All  the  / 

i= 1 

integrated  scores,  S' . ,  are  fused  (with  a  fusion  rule  that  can  take  the  form  of  a  logical  AND  or  OR 

operators,  or  a  majority  rule)  to  produce  the  system  declaration  “target  is  detected”.  A  system 
with  J  local  detectors  is  also  known  as  “distributed  detectors”,  “distributed  binary  integration 
system”,  and  “distributed  detection  network”.  The  challenge  is  to  properly  set  the  local  detectors’ 
double  thresholds. 

The  solution  for  combining  /  local  detectors  (each  with  N  scores)  as  a  distributed  binary 
integration  system  (also  known  as  a  double-threshold  detector,  or  m-out-of-n  detector)  is  given 
by  Han  et  al.10  where  the  cost  function  for  the  system  (i.e.  for  the  entire  car)  is  chosen  to 
minimize  the  system  probability  of  a  missed  detection,  Pmiss(system)  =  1  -  JPdetect(systeml ,  for  a  given 

constraint  of  constant  false  alarm  (CFAR)  system-false-alann  probability,  Pfalse(system)  ■  A  specific 

fusion  rule  (logical  AND  or  OR)  is  chosen  to  combine  all  the  J  scores  (each  score  is  produced 
with  a  binary  integration,  double  threshold  method).  The  method  of  solution  for  the  minimization 
problem  is  the  Lagrange  multiplier  method.  Let’s  assume  that  for  each  of  the /local  detectors  we 
have  N measurements  according  to  (2),  Mj=UoN  and  by  assuming  a  statistical  model  for  the 

signals  (nonnal  statistics  in  (10)),  one  can  compute  2  x/  detection  and  false  alann  probabilities 

oo  oo 

(from  (12))  that  apply  to  score n :  Pd(j)(x  >  //,  / )  -  J pl(x)dx  and  Plf]) (x  >  ijXj)  -  J  p(l (x)dx .  With 

hj  iij 

these  2  xj  quantitie  plus  the  system  false  alarm  probability,  Pfalse(smem),  and  /  x  N  scores 
(computed  from  Mi=UoN  for  the  /  local  detectors)  the  /  1st  threshold  //, .  values  and  /2nd 
threshold  rj2j  values  are  solved  for.  The  solution  process  involves  solving  simultaneously  2  x  / 

nonlinear  equations.  A  solution  is  obtained  with  numerical  optimization  routines. 

We  think  of  a  fusion  rule  in  the  fonn  of  logical  “OR”  as  the  proper  choice  for  detecting  a 
contaminated  car  that  is  observed  in  /  AORs  because  it  is  unlikely  that  all  the  AORs  are 
contaminated  (hence  we  rule  out  the  use  of  an  “AND”  operation)  and  also  a  majority  fusion  rule 
is  risky  because  it  is  possible  that  the  car  will  be  contaminated  only  in  one  AOR.  With  an  “OR” 
fusion  rule  we  declare  a  detection  if  at  least  one  (or  more)  of  the  /  local  detectors  (operating  with 
the  two  thresholds  binary  integration  method)  produces  a  detection.  Thus,  the  “OR”  fusion  rule 
is  more  likely  to  produce  better  system  detection  performance. 

A  full  implementation  of  the  Han  et  al.10  solution  for  distributed  sensors  is  not  easy 
(solving  2  x  /  nonlinear  equations),  and  thus,  we  give  an  alternative  solution  that  though  is  not 
as  optimal  is  much  more  practical  to  use.  Our  objective,  as  in  Han  et  al.  (1993)  is  to  set  the  two 
thresholds  r/Xj ,  rj2j  (for  each  of  the  /  local  detectors)  such  that  the  global  (system)  false 

alarm  Pfalseisystem)  is  set  to  be  CFAR  mode  of  operation,  and  to  compute  the  resulting  global 

detection  probability  Rdetect(iVStem) .  We  limit  ourselves  to  the  “OR”  fusion  rule.  For  simplicity,  we 
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also  set  all  the  J  thresholds,  ij2j ,  such  that  the  false  alarm  rate  for  all  j  are  identical.  In  (13-19) 
we  outline  or  solution  process  for  a  given  Pfahe(system)  ■  The  global  false  alann  for  an  OR  fusion  rule 
is  given  by 

j 

Pfalse(svstem)  ^  j  j  ^  P *  false,  7)  0^) 
j= 1 

where  Pfalse  j  is  the  false  alarm  of  the  jth  local  detector  (after  applying  the  second  threshold,  i.e.,  at 
the  second  stage  of  the  “double-threshold  detector”  method).  In  (13)  l-P/a/  .  is  the  probability 

j 

of  a  correct  H0  event,  and  the  product  j^[  (1  -  Pfalse  j )  is  the  joint  event  of  all  J  detectors  to  get  a 

j= i 

true  (correct)  Ho  event.  False  alann  is  the  complement  of  the  probability  of  correct  Ho  event,  i.e., 

1  -  (correct  H0) .  Each  of  the  false  alann  probabilities,  Pfalsej ,  is  given  by  the  binomial  distribution 
for  the  binary  integration  method  as 


N\ 


~ false,  j 


P^i\-PfljY-m  (14) 


m\(N-m)\ 

where  Pflj  is  the  false  alarm  of  the  /h  detector  operating  with  the  1st  threshold  (i.e.,  at  the  1st 
stage  of  the  double-threshold  procedure).  We  note  that  in  (14)  Pfahe  j  is  a  function  of  the  2nd 
threshold  tj2j  .  The  objective  is  to  solve  (13,  14)  for  all  J  lststage  Pflj  probabilities.  For 
simplicity  (and  hence  lack  of  optimality  in  the  system  performance)  we  assume  in  (13)  that 
all  Pfalse  j  are  equal,  hence, 


P  =  1  -  (1  -  P 

false  ,j  V  false(system) 


(15) 


We  substitute  (15)  in  (14)  and  solve  (with  a  numerical  solver)  for  Pf]/  as  a  function  of  the  2nd 
threshold  rj2j  .  Due  to  the  large  value  of  N  (N  can  be  100  or  more),  we  may  approximate  (14) 
with  a  Poisson  distribution  which  is  a  good  approximation  of  a  binomial  law  when  N  is  large  and 
Pflj  is  small.  Because  Pfahej  is  a  function  of  both  Pf]J  and  rj2j ,  we  need  to  choose  a  value  for 
ri2j  in  order  to  solve  for  Pflj  .  Currently,  we  pick  the  value  rj2J  that  maximizes  the  detection 
probability  Pd(J)  for  each  j  (this  improves  the  optimality  of  our  procedure,  but  still  does  not 
guaranty  an  optimal  overall  system  probability  of  detection  Edetect(  svstem) ,  hence  our  result  is 
suboptimal  with  respect  to  Han  et  ah).  We  use  (12)  to  solve  PPJ  for  the  1st  threshold  rjXj\ 


Pfij  =  J  P0  (x)dx  =  2~' [1  -  erf  ( 


2(*/l  tjSj) 


0] 


Vij  =  erf  'W  ~  2Pfij]j2(SjT  Zqj  8j) 


(16) 
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where  Z0/  and  8 .  are  the  covariance  (see  (10))  and  the  target  signature  (see  (11)),  respectively, 
for  the  jth  AOR.  Given  //l  ;  we  compute  the  1st  stage  local  detection  probability  with  (12)  as 


p«j  = 


ou 

j  /;,  (x)dx  =  2“'  [1  -  erf  ( 


Vi 


-1 


So  jSj) 


-2 


)]•  (IV) 


The  second  stage  local  detection  probability  is  given  in  a  similar  manner  as  (14)  with  the 
binomial  distribution  as  a  function  of  the  2nd  threshold,  ?/2  / ,  by 


^det ectj=  X  — ^-[PdlJP(Hl)ni-PdljP(Hl)]N'm  (18) 

where  P(Hf  is  the  prior  probability  for  ll\  event  (target  present  within  the  FOV).  When  the 
FOV  is  small  compared  to  the  footprints  of  the  contamination  (given  by  the  cross-sectional  area 
distribution  (pa(a)  in  section  3.3),  P{HX )  is  approximately  the  fraction  of  the  total  cross-sectional 

area  within  the yth  area  ( A/J in  (5)).  Exact  values  for  P(IIt )  can  be  obtained  via  the  simulation  for 


the  fill  factor  PAf )  discussed  in  section  3.3  and  Appendix  A.  In  general  we  seek  the  largest 


possible  threshold  that  still  produces  acceptable  probability  of  detection  due  to  the  fact  that  large 
threshold  improves  the  rejection  of  exceptionally  large  clutter  (unexpected  clutter  that  is  not 
predicted  by  the  model).  With  alU  local  detection  probabilities,  the  global  (system)  probability 
of  detection  for  an  “OR”  fusion  rule  is  given  in  a  similar  manner  as  (13)  by 

j 


detect  (system) 


^-fR-'W  as) 


j= i 


where  1  -  PdclcclJ  is  the  probability  of  a  missed  detection  and  the  product  (1  PdcKctJ )  is  the 

7=1 

joint  event  of  all  J  detectors  to  miss  a  detection.  Detection  is  the  complement  of  the  probability 

j 

of  a  miss,  hence  the  system  detection  probability  is  1  - JJ  (1  -  PdctcctJ  )  .  As  noted,  our  solution 

7=1 

process  is  practical  for  implementation  but  is  not  optimal.  We  can  attempt  to  optimize  our  choice 
for  r/2  j  by  repeating  (14-19)  for  different  values  for  rj2J  and  chose  the  value  that 

maximizes  Pdetect(system) . 

At  present  our  optimization  process  for  finding  the  two  thresholds  //, .  and  //2/  is  as 
follows.  Given  Pjalse(systen)  and  Jin  (15)  we  solve  for  the  false  alarm  Pfa,seJ  of  the  y'th  local 
detector.  We  use  Pfahej  in  (14)  and  find  Pflj  for  a  given  set-value  of  r/2  .  in  (14),  say,  //2/  =  1 . 
With  that  Pflj  we  can  compute  the  1st  threshold  rjXj  for  the /h  local  detector  in  (16),  and  we  also 
compute  the  1st  stage  detection  probability/^^  in  (17).  We  repeat  the  computation  of  Pdlj  for  all 
values  r/2j  =  1  to  N  ,  and  choose  the  specific  value  rj2J  that  maximizes  Pdlj  (this  is  the  best  value 
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for  the  2nd  threshold).  For  this  “best  2nd  threshold”  we  already  know  the  value  of PfXj ,  hence  we 
know  the  “best  1st  threshold”  77 1;  given  by  (16).  While  not  100%  optimal  with  respect 
to  detect ,  due  to  the  constraint  of  same  Pfake  j  for  all  j  in  (15),  our  selection  process  for  the 
two  thresholds  works  well  in  our  simulations. 

7.  System  Requirements 

We  address  the  scenario  of  a  contaminated  car  with  fill  factor/ and  contamination  thickness  h 
(both/ and  h  are  stochastic  variables  that  have  fluctuations  given  by  pdfs  in  section  3),  subject  to 
the  geometry  of  the  problem  (see  (2))  and  the  external  constraints  and  sensor  parameters  (see 
(5)).  The  geometry  parameters  are:  the  distance  of  the  laser  source  to  the  car,  rsource ;  the  distance 
from  the  contaminated  car  (the  target)  to  the  detector,  rdeteclor ;  and  the  atmospheric  transmission, 
talm  ,  that  is  a  function  of  range  (rsource  +  rdetector) .  Sensor  parameters  are:  total  time  ( T)  available 
for  scanning  a  suspected  car;  sensor  acquisition  rate,  PRF  (i.e.,  measuring  p  spectral  bands 
every  PRF1  seconds);  sensor  field  of  view,  FOV;  and  the  J  areas  of  regard  with  associated 

j 

viewed  area  A  =  y A ,  .  The  average  surface  density  for  the  contamination  G(g/m2)  is  given  by 

M 

the  product  (multiplication)  of  the  contamination  intrinsic  material  density  (2.34  g/cm  for 
potassium  chlorate)  and  the  average  thickness  of  h,  E(h).  This  gives  the  average  surface-density 
of  the  contamination  that  is  encountered  within  the  FOV  under  the  If  scenario.  The  global 
contamination  surface-density  within  the  viewing  area  Aj  is  much  smaller  (the  contamination  is 
not  spread  unifonnly  over  the  car’s  surface)  and  is  given  by  the  product  of  G(g/m  )  and  the  prior 
probability  P{H\).  We  compute  the  number  of  areas  of  regards  J,  number  of  samples  N,  and  the 
measurements  Mi=UoN  (see  (2)).  The  numerical  value  of  J  is  detennined  either  by  subjective 

requirement  such  as  “need  to  sample  the  car  at /specific  locations”,  e.g.,  J=  5  for:  hood,  two 
front  doors  and  trunk  area;  or  by  the  constraint/  =  A  /( N  x  r\ource FO V ) ,  where  N  is  the  number  of 

measurements  (e.g.,  number  of  pixels  in  a  2D  imaging  sensor)  and  FOV  is  the  instantaneous  field 
of  view  of  a  pixel.  With  the  solution  (section  6)  subject  to  the  system’s  requirements  of  a  given 
value  of  CFAR,  we  obtain  the  two  thresholds,  and  we  use  the  1st  threshold  .  as  follows.  The 

numerical  value  of  the  two  thresholds  includes  the  effect  of  clutter  via  the  uncertainty  that  is 
captured  in  the  statistics  of  Rq.  In  most  detection  scenarios  clutter  noise  is  much  larger  than  the 
sensor  noise  (noise-equivalent  spectral  radiance,  NESR),  we  expect  that  state  of  the  art  sensors 
(NESR  ~10"  W/cm  Isrlcm  )  is  sufficient.  A  system  design  involves  a  tradeoff  between  many 
parameters  (data  acquisition  parameters  and  scenario  parameters).  Our  stochastic  model  for 
contaminated  surfaces  can  help  in  studying  these  tradeoffs.  For  example,  given  a  requirement  of 
CFAR  operation  we  can  use  our  model  to  compute  detection  thresholds  (?ilj,?j2j)  as  a  function 

of  the  contamination  amount  ( h )  and  to  tie  the  numerical  value  of  the  thresholds  to  the  system’s 
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characteristic  noise.  In  addition,  with  our  stochastic  model  we  can  study  the  effect  of  stochastic 
fluctuations  in  various  parameters  (h,  f,  Ro)  on  the  detection  performance,  as  well  as  studying  the 
interaction  between  the  external  system  parameters  ( T,  PRF,  FOV)  and  the  detection 
performance. 

We  view  the  utility  of  our  stochastic  model  for  rough  surface  contamination  in  two  types 
of  scenarios  with  regard  to  the  evaluation  and  the  design  of  a  detection  system.  In  the  first 
scenario  we  evaluate  the  perfonnance  of  a  given  sensor  in  a  detection  scenario  where  the 
contamination  parameters  are  unknown.  Thus,  in  this  scenario  the  external  constraints  (5)  are 
given  (i.e.,  a  sensor  with  FOV,  and  PRF)  and  the  contamination  variable  ( h,f  Ro)  in  (2)  are 
variables.  Please  note  that  in  (5)  the  detection  is  only  based  on  the  total  number  of  measurements 
J  x  N  (where  each  of  the  N  measurements  contain  p  spectral  bands)  and  that  N  can  be  achieved 
by  measuring  the /h  AOR  with  an  imaging  system  that  collects  N  pixels,  or  a  sensor  with  a  single 
pixel  that  collect  N  measurements.  Thus,  if  the  FOV  of  a  system  and  the  range  to  the  suspected 
car  is  such  that  for  a  given  j  (c.g.,/=  I  is  the  front  door  and  the  sensor  is  at  a  distance  r  from  the 
car)  is  such  that  the  illuminated  footprint  FOV  x  r 2  on  the  car  is  equals  to  the  door’s  area,  then  N 
can  be  achieved  by  sequential  measurements.  On  the  other  hand,  if  the  FOV  is  such  that  there  are 
N pixels  within  the  door’s  area,  then  with  an  imaging  system  “one”  measurement  can  provide  the 
N  samples  in  (2,  5).  We  remind  the  reader  that  the  fill  factor  is  a  function  of  the  sensor’s  FOV  as 
well  as  a  function  of  the  contamination  cross-sectional  foot-print  and  the  contamination  level. 
With  our  model  we  can  estimate  the  sensor  detection  probability  as  a  function  of  the 
contamination  level.  The  sensor  intrinsic  external  constraints  {FOV,  PRF)  together  with  the 
constraints  about  the  detection  scenario  ( T,  A,  J)  are  assumed  to  be  given.  With  our  model  we 
can  estimate  the  detection  probability  Pdetect(h,  f ,  R0  \  CFAR,FOV,PRF,T,  A,J)  as  a  function  of 

the  contamination  thickness  (h),  the  fill  factor  (j)  and  the  rough  surface  reflectivity  (Ro),  for  a 
given  false  alann  (CFAR)  and  the  external  constraints  (FOV,  PRF,  T,  A,  J).  The  second  scenario 
for  the  utility  of  the  model  is  to  answer  the  following  question:  “What  should  be  the  best  sensor 
parameters  (FOV,  T,  PRF,  A)  for  given  a  contamination  scenario  with  (h,f  Ro)?”  For  this  utility 
scenario  we  can  compute  the  detection  probability  PAetect(FOV ,PRF ,T ,A,J  \h,f,R0).  Setting  a 

value  for  CFAR  depends  on  operational  considerations.  For  example,  if  the  scenario  is  of  a  car- 
stop  checkpoint  where  a  car  is  inspected  with  a  sensor  and  when  an  initial  detection  is  declared, 
the  car  is  diverted  for  a  further  (more  thorough)  check.  If  the  team  at  secondary  checkpoint  can 
handle  a  load  capacity  of  checking  one  car  out  of  20  (5%),  then  the  CFAR  for  the  initial 
checkpoint  (14-19)  can  be  set  to  CFAR= 0.05.  The  larger  the  CFAR  is,  the  higher  the  detection 
probability  is,  and  thus  the  probability  of  missing  the  presence  of  a  contaminated  car  is  reduced 
(which  is  the  primary  goal  of  the  checkpoint). 

There  are  three  key  parameters  that  most  effect  the  overall  detection:  the  number  of 
measurements  (N),  the  prior  probability  P(H\)  of  viewing  a  contaminated  area,  and  the  pseudo 
SNR  parameter  y  (given  in  (12))  which  is  related  to  the  strength  of  the  U\  scores  (i.e., 

/increases  with  increased  target  signal  and  is  inversely  proportional  to  the  variance  of  the  rough 
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surface  X0 ).  The  lower  P(H\)  is,  the  larger  A  should  be  in  order  to  enhance  the  likelihood  of 
intercepting  the  contamination  present  in  A  with  some  of  the  laser  shots.  The  larger  y  is,  the 
higher  the  detection  probability  is.  We  explore  the  interaction  between  the  different  variables  in 
section  8. 

8.  Results 

We  explore  the  parameter  space  ( h,f  Ro,  A,  J,  P(H\ ))  with  regard  to  probability  of  detection, 
presented  in  section  8.1,  and  also  present  a  case  studies  giving  the  system  detection  perfonnance 
for  detecting  specific  contaminated  cars  (section  8.2).  In  all  our  simulations  a  zero-mean 
additive  Gaussian  noise  with  standard  deviation  of  0.005  was  added  to  the  simulated 
measurements  in  (2).  The  objective  system  false  alann  was  set  in  all  simulations  to  be 
CFAR=0.05.  The  geometry  for  the  simulations  assumes  backscattering  angle  9  =  0/  in  (1),  and 

Ksystem=l  in  (2).  We  note  that  most  of  the  fluctuations  in  the  simulated  measurements  are  due  to 
the  stochastic  nature  of  the  RVs.  Additive  noise  is  usually  easy  to  handled  in  detection 
algorithms  because  it  is  white  (spectrally  isotopic).  The  challange  in  detection  algorithm  is  how 
to  handle  “structured”  noise  in  the  form  of  clutter,  which  is  given  in  our  work  by  the  pdf  of  Ro. 
Our  contamination  is  a  potassium  chlorate  material  (specific  density  of  2.34  g Icm  ).  In  all 
simulations  gAn  in  (1)  is  set  to  zero  for  simplicity,  the  contamination  is  potassium  chlorate,  and 
the  extinction  coefficient  a  was  taken  from  the  KBr  pellet  measurements  in  reference  1 .  Note 
that  the  RV/(fill  factor)  is  conditioned  to  the  II \  event,  i.e.,  it  is  the  fraction  of  contaminated 
area  within  the  illuminated  spot  (we  assume  that  the  FOV matches  the  laser’s  beam  divergence) 
when  contamination  is  actually  present.  The  prior  probability  for  H\  event  is  (roughly)  the 
portion  of  the  measured  area  that  is  contaminated,  where  the  measured  area  includes  subareas 
that  are  not  contaminated.  The  parameter  space  is  >  3  and  thus  it  is  not  possible  to 
simultaneously  display  the  intricate  interaction  between  the  variables  (for  this  we  need  a  “7D 
plot” — 6  independent  parameters  plus  a  probability  of  detection  as  a  dependent  variable).  We 
remind  the  reader  that  with  the  external  constraints  in  (5),  one  can  relate  physical  parameters  of 
the  sensor  (FOV,  PRF,  and  distance  to  the  surface),  and  acquisition  parameters  (total  area,  A,  to 
sample  and  time,  T,  to  acquire  the  measurements)  with  the  detection  parameters  (A  and  J).  For 
example,  given  A=100  and  J=4,  e.g.,  for  two  doors  (j=  1,  2),  hood  and  trunk  (/= 3,  4);  or  j=  1,  2  for 
painted/unpainted  vertical  surfaces,  and  j= 3,  4  for  painted/unpainted  horizontal)  and  a  sensor 
with  1  mr  linear  field  of  view  (solid  angle  FOV  is  10'6  sr)  at  a  distance  of  20  m  (i.e.,  2  cm 
illuminated  spot)  and  PRF  of  10Hz,  implies  that  total  time  for  measuring  the  contaminated  car  is 
T  =  40.v  =  J  x  N / PRF ,  and  the  total  car’s  area  that  is  measured 

is  ^4  =  1600  cm2  =  Ax  r2FOV  x  J  (i.e.,  a  small  10  cm  by  10  cm  area  for  each  of  the  4  J’s). 
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8.1  Exploring  the  parameter  space 

We  explore  the  parameter  space  (Ro,  h,f)  as  is  used  in  our  2nd  order  statistics  (mean  and 
variance  in  equations  9-12)  and  its  effect  on  detection.  We  compute  (13-19)  the  double 
thresholds,  the  overall  global  (system)  detection f,detection(ca/^) .  We  present  a  few  3D  plots  to  show 
the  complex  interaction  between  the  key  parameters  that  govern  the  detection  scenario: 

(^detect(system)  > X  N) .  We  constrain  the  desired  overall  (global)  system  CFAR  to 

be^faiSe-aiarm(car)  =  0-05  .  To  explore  the  parameter  space  we  selected  different  combinations  of 
mean  values  for / and  h,  E if)  and  E (h),  respectively,  and  arbitrarily  set  the  variances  as  follows: 
the  standard  deviation  of  h  is  10%  of  its  mean  value;  the  standard  deviation  (std)  of  /is  the 
minimum  value  between  E(f)/5  and  [  1  -E(f)]/5  (where  the  value  5  is  5-standard  deviations)  to 
ensure  that  for  normally  distributed  /  the  constraint  0  <  /  <  1  is  met.  Our  choice  for  standard 
deviation  of /  enforces  a  low  standard  deviation  near  the  endpoints  (f=  0  and/=  1 )  and  a  high 
standard  deviation  in  the  middle  (/=0.5).  Our  choice  for  std  can  be  interpreted  as  follows:  a  low 
std  near /=0  may  occur  for  a  sensor  with  large  FOV viewing  sparse  contamination  with  a  small 
cross-sectional  area;  a  small  std  at/=l  may  occur  for  a  sensor  with  a  small  FOV  and 
contamination  with  large  cross-sectional  area;  and  for  y~0.5  the  sensor  FOV  matches  the  size  of 
the  contamination  footprint  and  how  the  FOV  and  contamination  overlap  (intersect)  can  produce 
both  large  and  small  values  of  /  hence  a  large  std.  The  mean  of  Ro  is  taken  from  for  lab 
measurements1  of  rough  aluminum  where  a0  =  6.96  x  1 0  8  cm 2 ,  b0  =  -6.98  x  10  4  cm  ,  c0  =  0.25 

and  the  corresponding  standard  deviations  were  arbitrarily  set  at  (lxl 0  7 cm2 .1x10  8 cm,  lxl  0  3 ) 
as  we  do  not  have  enough  lab  data  to  estimate  standard  deviations.  Elpon  availability  of  data  we 
can  update  the  standard  deviation  values.  The  Fresnel  coefficient  p0  is  for  pure  aluminum 

(whose  reflectivity  across  the  LWIR  is  -0.99).  Due  to  the  large  aluminum  reference  the  effect  of 
surface  reflection  from  the  target  was  small  and  allowed  us  to  set  a,  =b,  =  c,  =  - oo  in  (1).  The  19 
(p)  wavenumbers  in  our  simulations  correspond  to  available  CO2  laser  wavelengths  (9-1  lp m), 
given  in  cm'1  as  934.93,  942.42,  949.49,  956.21,  969.18,  974.66,  975.9,  978.47,  982.13, 
986.58,  1031.5,  1039.4,  1046.9,  1052.2,  1057.3,  1070.4,  1074.7,  1078.6  1083.4. 

In  Fig.  la  we  show  the  interaction  between  ( /ietcct(System) ,y,P(Ht),  N)  for  a  given 

CFAR=0.05  and  for  one  region  (J=  1 )  of  aluminum  rough  surface.  For  each  combination  of  E(/) 
and  E (h)  the  pseudo-SNR  parameter  y  was  computed,  and  a  series  of  system  detection 
probabilities  were  computed  for  different  values  of  N  and  P (//).  The  figure  shows  that  the 
overall  system  detection  probability  increases  as  P{H\),  N,  and  y  increase.  In  Fig.  lb  we  show 
the  relationship  between  y  and  the  means  of  the  fill  factor/ and  the  contamination  thickness  h 
(may  be  converted  to  surface  contamination  density  via  G(pg  /  cm2)  =  234 x  E(h)  where  h  is  in 
microns).  The  figure  shows  the  mono  tonic  increase  of  y  with  h  and  / 

We  can  also  use  the  stochastic  model  to  set  sensor  requirements  for  a  contamination 
scenario.  In  Figs.  2a  and  2b  we  show  an  expanded  portion  of  Fig  1.  Note  that  Fig.  2a,  for 
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example,  can  be  used  to  determine  what  values  for  parameters  N,  P{H\),  and  y  are  required  to 
achieve  a  particular  system-level  probability  of  detection.  For  example,  assuming  a  desired 
P detect( system)  >  0.8,  any  point  lying  on  or  above  the  red  contour  surface  will  satisfy  the  desired 
condition  (given  the  pdf  choices  used  to  create  the  figure).  The  mapping  is  not  unique  (many 
sets  of  values  for  N,  P{H\)  and  y  will  satisfy  the  condition)  and  therefore  it  is  difficult  to  invert 
Eqs.  (13-19)  to  solve  for  N,  P(JI\)  and  y  given  /Ticiccl(-syslcm).  However,  desired  information  about 

N,  P{H\),  and  y  given  Pdetect(system)  can  be  read  off  of  plots  such  as  Fig.  2a.  The  parameters  N, 
P{H\),  and  y  all  relate  to  sensor  requirements  (external  constraints  in  Eq.  (5)) — -N  relates  to  the 
time  to  detect  and  FOV  of  the  sensor,  P(H\)  depends  on  the  FOV  and  the  nature  of  the 
contamination,  and  y  depends  both  on  the  nature  of  the  contamination  (thickness  h,  fill  factor/, 
surface  density  G)  and  the  sensitivity  of  the  sensor  (laser  power,  noise  level,  etc.). 
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Fig.  la  System  detection  probability  for  one  area  of  regard  (2=1)  of  rough  aluminum  as  a  function  of  number  of 
measurements  (A),  prior  probability  P{H\),  and  the  pseudo  SNR  parameter  v .  Contour  surfaces  show  where 

P detect! system)  =  constant.  Three  contour  surfaces  are  shown  at  Pdetect(system)  values  of  0.2  (blue),  0.5  (green/yellow),  and 

O. 8  (red). 


E(h)  [pm]  E(f) 

Fig.  lb  Pseudo-SM?  parameter  y  as  a  function  of  the  mean  fill  factor,  E(f),  and  the  mean  contamination  thickness, 
E(h),  where  the  surface  contamination  density  is  G( pg  /  cm  ~  )  =  234  x  h  where  h  is  in  microns). 
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Fig.  2b  An  expanded  view  of  Fig.  lb  to  show  details  of  pseudo  SNR  parameter  y  as  a  function  of  fill  factor  and 
contamination  thickness.  Combinations  of  E(/?)  and  E (/)  resulting  in  y> 4  are  shown  as  dark  red  pixels,  y< 4  are 
mapped  to  other  colors  as  shown. 

In  Fig.  3  we  show  the  improved  detection  when  3  rough  aluminum  AOR  (J=  3)  are  used 
for  detection  (19)  instead  of  only  one  AOR.  The  improved  detection  can  be  observed  by 
comparing  Fig  3  with  Fig.  2a  where  it  is  evident  that  the  contours  for  high  detection  probability 
(brown  and  red  colors)  in  Fig.  3  spread  to  the  direction  of  lower  values  of  N,  P(H\),  and  y . 


Fig.  2a  An  expanded  view  of  Fig.  la  to  show  the  behavior  for  y< 4. 
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Fig.  3  Same  as  Fig.  la  but  J=  3  areas  of  regard  are  used  for  detection  in  (19). 

We  explore  the  effect  of  a  different  rough  surface  material,  painted  aluminum,  which 
may  be  more  realistic  for  a  car.  We  seta0  =  -3.47  x \0  7 cm2 ,  b0  =  -4.78  x  lCT4cm  ,  c0=1.07 
(derived  from  lab  measurements1)  and  the  corresponding  standard  deviations  are  the  same  as 
those  previously  selected  for  rough  aluminum  (1x10  7cm2, 1x10  5cm,  1x10  3).  Due  to  the  low 
reflectivity  of  the  painted  aluminum  ( p0  ~  0.07  across  the  19  laser  wavelengths),  the  potassium 
chlorate  rough-surface  parameters  (see  equation  1)  were  used  and  set  to  at  =  -3.20  x  10  1  cm2 , 
bt  =  2.07  x  1 0  3 cm  ,  ct  =  -2.58  (derived  from  a  lab  measurements1).  Corresponding  standard 
deviations  were  set  arbitrarily  to  (5  xlO  9cm2,5x  \0  5cm,5x  1 0  3 ) ,  as  we  do  not  have  enough  lab 
data  to  estimate  standard  deviations. 

We  observed  that  the  painted  aluminum  surface  produced  less  signal  (Min  equation  1) 
and  that  the  pseudo  SNR  parameter  y  was  about  factor  3  reduced  from  the  rough  aluminum 
surface  (in  Fig.  lb).  Thus,  the  system  detection  probability  will  be  reduced  given  the  same  E(h) 
and  E(f)  values.  One  may  use  Figs  la  (or  2a)  to  obtain  a  value  of  interest  for  ;/  that  produce  the 
desired  system  detection  probability,  and  find  in  Fig.  lb  (or  2b)  the  values  for  E(h)  and  E(f)  that 
corresponds  to  a  value  of  3 y  .  These  values  of  E(h)  and  E(f)  are  the  one  that  corresponds  to 
contaminated  painted  aluminum  surface. 

8.2.  A  case  study:  the  detection  of  a  contaminated  car 

We  explore  the  detection  scenario  with  two  case  studies  (simulation  1  in  Figs.  4-6,  and 
simulation  2  in  Figs  6-7)  for  a  contaminated  car  at  a  distance  of  20 m  to  be  detected  in  10s  where 
two  areas  of  regard  (J=  2)  are  sampled:  the  first  AOR  is  rough  aluminum,  and  the  2nd  AOR  is  a 
painted  rough  aluminum  (the  same  parameters  for  a0,  b0,  c0  presented  in  Section  8.1  were  used 
here).  The  pdf  for  h  is  a  ramp  (consistent  with  our  expectation  that  a  high  target  transmission  is 

h  -h 

more  likely  than  a  low  transmission)  and  thus,  Ph  (h)  =  2 - — - -  for  hmm  =  0.75 pm  and 

(/2,nax  ~K\n) 
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hm  in  =  0.0  lion  (Section  3.2).  For  this  distribution  E(h)  =  0.25  pm  .  For  the  pdf  for  /  we  present 
two  case  studies.  First  case  is  computed  (Appendix  A)  from  cross-sectional  area  (footprint  or 
spill  distribution),  <pa ,  measured  by  MIT-LL  (large  majority  of  the  footprints  are  small  for  which 

the  key  parameter  P(H\)  is  large — about  0.4 — see  pdf  of  fill  factor  and  P{H\)  in  Fig  Al).  For 
case  2  we  simulate  the  contamination  footprint  (Fig.  A2)  such  that  large  footprint  are  more  likely 
in  cpa  for  which  the  key  parameter  P(II\ )  is  small  (0.07).  In  our  scenario  the  laser  linear  FOV 


(and  beam  divergence)  is  1  mr  (i.e.,  a  2 cm  spot  size  on  the  car)  and  the  laser  is  configured  in  the 
traditional  backscatter  lidar  configuration.  The  system  constant 


K 


system 


source  (  / ) 


(fi 


4-  y 

source  detector ) 


is  set  to  be  one  in  (2).  The  laser  system  operates  at  10  Hz  and 


acquires  N=  50  spectral  measurements  for  each  of  the  AORs  (each  of  the  N  measurements 
contains  19  wavelengths  in  the  LWIR  range).  The  total  acquisition  time  for  J=2  is 
T  =  10  s  -  JxN  /  PRF  as  required. 

In  Fig.  4  we  show  results  for  case  1  for  the  effect  of  the  second  binary  threshold 
ij2j  (given  in  equation  18)  for  the  two  local  detectors  (j=  1  and  j= 2)  on  the  local  detection 


probability  Pdetect  .  and  also  on  the  system- level  detection  probability  Pdetect(system)  (see  (19)).  The 

pseudo-SNR  parameter  y  equals  2.7  and  1.2  for  j=  1  and  j— 2,  respectively,  meaning  that  it  is  an 
easier  detection  problem  to  detect  potassium  chlorate  signals  on  rough  aluminum  rather  than 
painted  aluminum.  As  a  result,  the  /=  1  region  shows  better  performance  than  j= 2  (maximum 
value  of  detect, >2  =  0.926  whereas  Pdctcct  /=l  exceeds  this  value  for  rj2l  <  8  ).  The  system-level 


performance  shows  a  larger  probability  of  detection  at  the  same  value  of  r/2j ,  or  alternatively  for 

the  same  probability  of  detection  a  larger  threshold  may  be  used.  We  desire  the  largest  possible 
threshold  that  still  produces  high  probability  of  detection  because  high  threshold  improves  the 
rejection  of  exceptionally  large  clutter  (unexpected  clutter  that  is  not  predicted  by  the  model). 

In  Fig.  5  we  show  the  overall  system  detection  probability  Pdetect(system)  in  (19)  for  case  1  as 


a  function  of  all  possible  combination  of  thresholds  ri2]  (given  in  18).  The  asymmetry  in  Fig.  5 
is  interesting.  It  shows  that  it  is  more  important  to  optimize  tj2l :  given  an  optimal  value  of  rj2l ,  it 
is  less  important  (Pdetect(system)  is  less  sensitive)  to  the  value  of  r/22 .  This  makes  sense  since  j=  1  is 
a  better  region  (higher  perfonnance)  than  j= 2. 
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Fig.  4  Detection  probability  of  local  detectors  Pdctccl  .  in  (19)  for  simulation  1  (footprint  contamination  Fig.  Al)  as 
a  function  of  the  2nd  binary  threshold  >)-, . .  j=  1  is  AOR  for  rough  aluminum,  and  j=  2  is  AOR  for  painted  aluminum 
surface.  P(H{)=§A.  The  system  detection  probability,  under  the  condition  that  r)2l  =  r/-n  ,  is  also  shown  in  black. 
PdctcctJ,2  reaches  a  maximum  of  0.926  at  ij22  =  2.  Red  dot  shows  the  largest  ij2l  such  that  Rdetcct;=1  exceeds  a 
detection  probability  of  0.999  (occurs  at  r/2]  =  3  ).  Black  dot  shows  the  location  where  the  system  detection 
probability  exceeds  0.999  (occurs  at  r/2l  =  tj22  =5).  Combining  the  two  regions  allows  larger  l]2j  thresholds  to 
achieve  the  same  detection  probability  and  will  result  in  more  robust  performance. 


Fig.  5  Overall  system  detection  probability  Pdetect(system)  in  (19)  for  case  1,  as  a  function  of  all  possible 
combinations  of  2nd  thresholds  772/  in(18).  P(Hi)=0A.  //,  is  2nd  threshold  for/=l  (rough  aluminum)  and  rj22  is  2nd 
threshold  for  j=2  (painted  aluminum).  Black  line  occurring  along  the  diagonal  (where  ij2l  =  //,, )  is  the  same  black 
line  appearing  in  Fig.  6;  red  and  green  lines  appearing  at  t]22  =  N  and  //2I  =  N ,  respectively,  are  the  same  as  in 
Fig.  4  (at  rj^j  —  N  the  detector  will  never  alarm  and  thus  information  from  the  jth  region  is  discarded).  Maximum 

of  ^detect(system)  occurs  at  /?21  =  1,  J]22  =  2  at  a  value  of  0.9999996. 
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In  Figs  6-7  we  present  the  detection  for  case  2  (Fig.  A2,  with  low  P(H\)=0.07).  The 
combination  of  thresholds  (/=  1  and  j= 2)  is  advantageous  to  produce  a  higher  overall  system 
detection  probability  compared  to  the  individual  regions.  Fig.  6  in  comparison  to  Fig.  4  shows 
the  importance  of  the  P(H\)  parameter,  where  in  spite  of  the  large  pseudo-STW?  parameter  y 
(equals  16.9  and  7. 1  for  j—  1  and  j= 2,  respectively)  that  is  much  larger  than  the  y’s  in  Fig.  4,  the 
overall  system  performance  is  worse  than  in  Fig.  4  because  P(H\ )  was  reduced  from  0.4  in  Fig.  4 
to  0.07  in  Fig.  6.  The  reduction  in  P{H\ )  means  that  the  contamination  is  more  sparse  and  target 
material  is  less  frequently  encountered  by  the  sensor.  The  parameter  P{H\)  plays  a  crucial  role  in 
system  performance  (i.e.,  detection  probability)  as  is  demonstrated  in  Figs.  1-3.  The  large  effect 
of  P{H\ )  can  also  be  seen  in  the  performance  of  the  two  individual  local  detectors  (/=  1  for  rough 
aluminum  and  j= 2  for  painted  aluminum).  In  Fig.  4  a  difference  of  1.5  units  in  y  (2.7  -  1.2) 
between  the  two  local  detectors  (red  and  green  lines),  results  in  large  difference  in  detection 
probability  (for  the  same  threshold),  while  a  difference  of  9.8  units  in  y  (16.9  -  7.1)  in  Fig.  6 
makes  very  little  difference  in  detection  probability  (the  red  and  green  lines  are  almost  identical). 

Note  that  in  the  second  simulation  the  best  performance  is  clearly  given  by  rj2l  =  rj21  =  1 , 
i.e.,  the  presence  of  a  single  non-zero  value  Sji  for  i  =  1,2  and  a  given  j  is  enough  to  declare 
that  the yth  region  is  contaminated.  This  does  not  mean  that  the  remaining  N~  1  zero-valued 
scores  were  “wasted”  and  that  one  could  have  set  N-  1  and  achieved  the  same  performance. 
Because  of  the  low  probability  of  F’(Hi),  a  large  N  is  required  in  order  to  increase  the  likelihood 
that  the  contamination  is  viewed  in  at  least  one  of  the  N  measurements. 


Fig.  6  Same  as  Fig  4  but  for  contamination  footprint  (case  2)  of  Fig.  A2.  />(//1)=0.07.  Individual  detection 
probability  Pdetect .  attain  a  maximum  of  0.934  and  0.919  for  j=l  and  j=2  respectively  at  rj2  ,  =  1 .  Maximum 

detection  probability  for  the  system  is  0.995  at  Tj2j  =  1  ■ 
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Fig.  7  Same  as  Fig.  5  but  for  contamination  footprint  (case  2)  of  Fig.  A2.  P(Hi)=0.01 . 

9.  Summary 

We  developed  a  framework  to  evaluate  the  performance  of  a  system  whose  objective  is  the 
detection  of  a  contaminated  car.  Our  approach  is  a  “system  approach”  where  external  constraints 
(5),  such  as  time  to  detect,  area  to  scan,  distance  to  the  car  are  input  parameters.  Our  algorithm 
employs  a  data  fusion  approach  known  as  a  distributed  binary  integration  system  (also  known  as 
a  double-threshold  detector,  or  m-out-of-n  detector)  in  order  to  combine  the  individual  detection 
results  from  multiple  scans  over  several  potentially  contaminated  areas.  Our  system  consists  of  J 
areas  of  regard  (e.g.,  hood,  doors,  etc)  that  are  scanned  with  an  active  laser  source  (p 
wavelengths).  Our  detection  algorithm  utilizes  the  matched  filter  concept  (which  is  optimal  for 
additive  Gaussian  noise;  in  our  detection  scenario  the  clutter  noise  is  not  additive)  together  with 
Neyman-Pearson  methodology  of  CFAR  operation  mode.  We  view  our  system  as  a  “distributed 
binary  integration  system”,  a  methodology  that  is  used  in  radar  literature.  Each  of  the  J  areas  is 
viewed  as  a  “local  detector”  and  the  system  fuse  the  local  detectors  output  (fusing  rules  can  be 
logical  AND,  OR,  or  majority  rule;  we  chose  the  OR  logic  in  this  report  in  order  not  to  miss  a 
detection)  into  a  “global”  detector  (i.e.,  “the  system”)  that  produces  a  declaration  “target  is 
detected”.  We  operate  each  of  the  detectors  with  the  “binary  integration,  double  threshold”  radar 
methodology,  a  methodology  that  is  robust  to  strong  clutter  fluctuations.  From  each  of  /areas  we 
obtain  N  spectral  measurements  (each  measurement  contains  p  spectral  bands).  For  a  given 
constraint  of  false  alarm  ( CFAR  mode  of  operation)  we  produce  the  desired  operational 
thresholds  for  each  of  the  J  local  detectors.  Our  solution  is  not  optimal  with  respect  to  overall 
detection  (or  missed  detection)  but  it  is  much  simpler  for  implementation  than  the  optimal 
solution  given  by  Han  et  al.10  where  he  minimizes  the  probability  of  a  missed  detection  for  a 
given  CFAR.  In  our  detection  scenario  the  car  is  viewed  as  a  rough  surface  that  is  partly 
contaminated  with  a  given  amount  of  contamination.  In  our  detection  model  all  the  scenario 
parameters  (contamination  thickness  h,  fill  factor f  and  the  car’s  reflectivity  Rq)  are  random 
variables  with  probability  functions. 
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We  think  that  building  a  detection  system  will  require  the  “snapshot  advantage 
architecture”  that  enables  the  acquisition  on  J  x  N  x  p  measurements  within  a  time  T  (traditional 
pushbroom/whiskbroom  types  of  acquisitions  are  much  less  efficient  and  produce  lower  signal- 
to-noise  ratios).  The  framework  that  we  developed  in  this  report  will  be  the  basis  for  developing 
a  system  approach  for  a  Raman  detection  system  and  for  a  passive  infrared  system  where  the 
incidence  on  the  contaminated  car  is  the  ambient  thermal  radiation.  A  system  design  involves  a 
tradeoff  between  many  parameters  (data  acquisition  parameters  and  scenario  parameters).  The 
stochastic  model  for  contaminated  rough  surface  with  external  system  parameters  is  a  useful 
model  for  designing  and  prediction  the  performance  of  a  detection  system.  Given  a  requirement 
of  CFAR  operation  we  are  able  with  our  model  to  compute  detection  thresholds  (7,  .,  72/  )  as  a 

function  of  the  contamination  amount  (h,  G ).  In  addition,  with  the  model  we  can  study  the  effect 
of  stochastic  fluctuations  in  various  parameters  ( h,f  Ro)  on  the  detection  performance,  as  well  as 
studying  the  interaction  between  the  external  system  acquisition  parameters  ( T ,  PRF,  FOV)  and 
the  detection  performance. 

With  our  stochastic  model  we  can  explore  the  probability  of  detection  for  a  given  sensor 
(e.g.,  specific  FOV,  PRF)  under  variable  contamination  scenario  (e.g.,  h,f ), 

^detect  (contamination  |  sensor) ,  and  the  probability  of  detection  for  a  given  contamination  scenario 
for  a  sensor  with  variable  properties  (e.g.,  FOV,  PRF,  T).  We  can  also  use  the  stochastic  model 
to  set  sensor  requirements  for  a  contamination  scenario.  For  example,  Fig.  2a  can  be  used  to 
determine  what  values  for  parameters  N,  P{H\),  and  y  are  required  to  achieve  a  particular  system- 
level  probability  of  detection,  say,  Rdetectcsystem)  >  0.8.  The  parameters  N,  P(H\),  and  y  all  relate  to 
sensor  requirements  (external  constraints  in  Eq.  (5))  —  N  relates  to  the  time  to  detect  and  FOV  of 
the  sensor,  P{H\)  depends  on  the  FOV  and  the  nature  of  the  contamination,  and  y  depends  both 
on  the  nature  of  the  contamination  (thickness  h,  fill  factor  f  surface  density  G  )  and  the 
sensitivity  of  the  sensor  (laser  power,  noise  level,  etc.). 

In  the  report  we  showed  an  example  of  exploring  the  parameter  space  (h,f  Ro,  N,  J, 

P(H\)  with  regard  to  probability  of  detection,  and  we  also  explored  the  pdf  space  for  the  different 
variables  and  recommended  what  pdfs  to  use  for  each  variable  for  physical  simulations.  The 
most  important  parameters  that  affect  the  overall  detection  probability  are  the  number  of 
measurements  (N),  the  prior  probability  of  P{H\)  of  the  contaminated  areas,  and  the  pseudo-SNR 
parameter  y  which  is  related  to  the  strength  of  the  II\  scores.  The  lower  P{H\)  is,  the  larger  N 
should  be  in  order  to  enhance  the  likelihood  of  intercepting  the  contamination  present  in  A  with 
some  of  the  laser  shots.  The  larger;/ is,  the  higher  the  detection  probability  is.  With  our 
probability  model  we  can  explore  the  parameter  space  and  study  the  tradeoffs  between 
parameters  that  affect  the  overall  system  detection.  We  presented  3D  plots  (section  8)  that 
demonstrate  the  interaction  between  parameters  and  an  example  for  the  detection  of  a 
contaminated  car  with  a  CO2  tunable  laser  (19  wavelengths)  system.  In  our  examples  a  clutter 
distribution  for  Ro  was  assumed.  Our  choice  may  not  accurately  represent  the  variability  of  the 
clutter  observed  in  practice.  Thus,  a  key  recommendation  is  to  better  characterize  both  the  clutter 
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and  target  variability  (for  example,  through  empirical  measurements)  so  that  our  performance 
estimates  may  be  updated. 

This  report  is  the  foundation  (framework)  on  which  we  will  develop  a  stochastic 
detection  model  for  passive  infrared  remote  sensing  system  and  for  active  Raman  detection 
system.  In  the  next  report  we  plan  to  do  the  following: 

•  For  Raman  detection  R{6, 0j )  in  ( 1 )  will  be  modified  to  comply  with  the  physics 
of  Raman  scattering. 

•  For  passive  remote  sensing  (2)  will  be  replaced  by  Eq.  (1 1)  of  [1]. 

•  We  will  complete  our  theory  for  the  probability  of  a  till  factor  based  on  the 
random  distribution  of  centers  in  a  hyper  volume,  and  the  geometry  of 
intersection  between  the  FOV shape  (rectangular,  circle,  and  ellipse  for  tilted 
surfaces)  and  the  contamination  footprint-shape. 

•  We  will  explore  the  effect  on  overall  system  perfonnance  (19)  by  using  a  target 

direction  8  — »  a  instead  of  8  — »  jU0  in  (1 1). 

•  We  will  optimize  our  solution  for  maximizing  system  performance  by  allowing 
each  of  the  local  detectors  in  (15)  to  have  separate  probability  of  false  alann. 

•  We  will  explore  the  quality  of  using  2nd  order  statistics  (mean  and  variance)  in 
(12)  with  histogram  computed  from  sampling  the  RVs. 

•  We  will  explore  the  probability  space  of  the  RVs  (i.e.,  the  parameters  that  control 
the  shape  of  the  pdj s  in  section  3)  and  its  effect  on  system  detection  probability. 

Appendix  A.  Fill  factor  probability 

We  estimate  the  fill  factor  of  the  contaminated  by  a  brute  force  approach  where  we  simulate  the 
random  contamination  on  the  surface  as  produced  by  a  given  size  distribution^  (a)  of  the 
contamination  events.  Within  each  spot  the  material  thickness  (h  in  equation  1)  is  not  necessarily 

3 

uniform,  but  the  mean  surface  density  of  the  material  is  given  by  G’=2.34xE(/?),  where  2.34  g/cm 
is  the  intrinsic  density  for  potassium  chlorate  E(h)  is  the  average  height  of  the  contamination 
(i.e.,  G  is  the  average  contamination  surface  computed  over  all  of  the  “spots”  in  Figs.  A1  and 
A2).  The  sparsity  of  the  contamination  “spots”  in  Figs.  A1  and  A2  is  given  by  the  parameter 
P(Hi )  in  (18)  is  a  key  parameter  that  affect  system  performance  (see  section  8).  Accurate 
computation  of  P(H\)  is  affected  by  the  size  of  the  sensor’s  FOV,  because  sparsity  is  “measured” 
with  regard  to  FOV  (details,  below).  Figs.  Aland  A2  are  used  in  section  8.2  as  case  studies  for 
contaminated  areas  of  a  car.  The  general  idea  behind  the  simulations  is  that  we  assume  that  the 
locations  of  the  target  particles  (their  centers)  are  uniformly  distributed  on  the  surface7.  Then,  we 
randomly  placed  a  square  field  of  view  with  area  A fov  (m  )  on  the  simulated  contaminated 
surface  and  count  the  area  At  (m  )  for  the  target  (contamination)  within  the  field  of  view.  The 
ratio  At  A4fov  is  the  fill  factor. 
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Details  of  the  simulation  are  as  follows.  Let  a  denote  the  cross-sectional  area  of  a  single 
contamination  event  (i.e.,  a  fingerprint,  a  single  particle,  a  thin  film,  or  a  spill  of  many  particles). 
We  make  the  following  assumptions:  a  surface  of  with  total  area,  A,  is  potentially  contaminated 
with  multiple  contamination  events.  The  contaminated  area  of  all  target  events  is  At,  meaning 
that  the  prior  probability  of  contamination  is  At/A;  the  likelihood  that  any  infinitesimal  portion  of 
the  total  area  is  covered  by  contamination  is  uniform;  the  likelihood  that  a  single  contamination 
event  is  of  a  certain  cross-sectional  area  is  given  by  pdf  Pa(a )  (computed  from  (pa(a)  in  section 

3.3),  contamination  events  are  assumed  to  be  circular.  Areas  from  multiple  contamination  events 
may  not  overlap  (reasonable  for  sparse  contamination). 

Under  these  conditions,  a  simulated  surface  is  populated  in  a  multi-step  process.  The 
simulated  surface  is  a  digital  image  where  each  pixel  in  the  image  is  of  area  Apixei  .The  digital 
image  is  initialized  with  each  pixel  taking  a  value  of  zero.  At  the  end  of  the  generating  process, 
the  digital  image  will  contain  a  value  of  the  area  within  each  pixel  that  is  covered  by 
contamination.  The  simulation  steps  are  as  follows: 

1 .  Determine  average  number  of  contamination  events  needed  to  achieve  At:  N,  =  E (a)/At 
where  E(a)  denotes  the  average  contamination  event  area  (expected  value). 

2.  Sample  N,  contamination  target  events  from  pdf  Pa  (a)  to  obtain  ah  i=l,2,...,Nt.  For 

convenience,  sort  the  contamination  events  in  order  of  descending  area. 

3.  Place  the  7th  contamination  event  by  choosing  a  pixel  at  random  which  will  serve  as  the 
center  of  the  contamination  event,  excluding  pixels  that  would  cause  overlap  between 
contamination  events.  For  the  first  contamination  event,  no  pixels  are  excluded.  Note 
that  contamination  events  smaller  than  a  single  pixel  may  be  placed  into  an  already 
contaminated  pixel  if  the  area  of  contamination  in  that  pixel  is  1-  ai/Apixei  or  less,  though 
this  will  happen  with  a  reduced  probability. 

4.  If  a,  >  Vati Apixei  the  contamination  event  must  span  multiple  pixels  (the  fraction  z/4 
accounts  for  the  difference  in  linear  dimension  between  a  circular  contamination  and  a 
square  pixel),  and  the  event  is  placed  in  the  image  using  a  Euclidean  distance  transfonn: 
all  pixels  within  a  distance  of  n  units  will  be  assigned  an  area  of  Apixe/,  where  r,  is  the 
corresponding  radius  in  pixels  of  the  contamination  event  if  it  were  circular.  Note  that 
this  step  is  an  approximation,  since  edge  pixels  around  the  contamination  would  more 
accurately  be  represented  with  values  less  than  Apixei. 

5.  If  a i  <  Vai iAPiXei  the  contamination  event  can  be  contained  entirely  within  a  single  pixel, 
and  a  value  of  ailApixei  is  added  to  the  current  value  of  the  pixel.  This  additive  adjustment 
is  the  increase  the  contaminated  area  within  the  pixel  due  to  the  additional  contamination. 
Note  that  this  step  ignores  the  fact  that  a  contamination  event  may  not  necessarily  occur 
in  the  exact  center  of  a  pixel  and  could  still  span  over  multiple  pixels. 

6.  Once  all  contamination  events  are  placed  in  the  image,  the  sensor  FOV is  translated 
across  the  image  computing  the  total  contaminated  area  within  the  FOV  at  each  FOV 
location;  the  fill  factor / as  seen  by  a  sensor  is  obtained  by  dividing  by  A Fov- 
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A  histogram  of/ gives  an  estimate  of  Pf{f)  as  seen  by  a  sensor;  i.e.,  it  gives  the  probabilities 

for  the  sensor  to  observe  a  specific  value  of  /  In  addition  the  simulation  can  give  P(H\ ) 
given  by  the  ratio  of  number  of  FOV s  that  are  not  empty  (i.e.,  FOV s  that  contain  any  target 
material)  to  the  number  of  FOV s  (A/Afov)  that  are  contained  in  the  area  A.  For  a  “small” 
FOV  (e.g.,  an  FOV  that  is  smaller  than  the  size  of  the  contamination  “spots”  in  Figs  A1  and 
A2),  P{H\)  is  simply  given  by  the  fraction  of  the  area  occupied  by  the  target  material,  A,  /A. 
Flence,  P(H i)  can  be  viewed  as  a  global  (macro)  property  of  the  contaminated  surface,  and  G 
is  the  local  (micro)  property  of  the  contamination  surface  density  of  the  contamination  events 
(“spots”). 


Fig.  Al.  Simulation  1  (case  1):  (left)  value  of/for  each  pixel  on  the  simulated  surface.  The  sensor  FOV  (blue 
box,  l/10th  the  linear  dimension  of  the  image)  moves  over  the  contaminated  surface,  (right)  The  fill  factor  as  seen 
by  the  FOV  for  the  contaminated  surface.  The  likelihood  for  the  sensor  FOV  to  encounter  contamination  is 
P(Hi)=  0.4. 


Fig.  A2  Simulation  2  (case  2):  same  as  Fig.  Al  except  that  individual  contamination  events  are  larger  in  cross- 
sectional  area;  it  takes  fewer  events  for  the  same  total  area  of  contamination.  The  likelihood  for  the  sensor  FOV  to 
encounter  contamination  is  P(//i)=0.072. 
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Fill  factor  (f) 


Figure  A3.  Fill  factor  probability  for  simulations  in  Fig  A1  and  A2. 
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